AD/A-004  479 

APPLICATION  OF  ADVANCED  METHODS  FOR  IDENTIFICATION 
AND  DETECTION  OF  NUCLEAR  EXPLOSIONS  FROM  THE  ASIAN 
CONTINENT 


T.  C.  B?che,  et  al 

Systems,  Science  and  Software 


Prepared  for: 

Air  Force  Office  of  Scientific  Research 
Advanced  Research  Projects  Agency 


October  1974 


DISTRIBUTED  BY: 


National  Technical  Information  Service 
U.  S.  DEPARTMENT  OF  COMMERCE 


UNCLASSIFIED  — 


SECURITY  CLASSIFICATION  CF  This  PACE  (Hhen  Dnla  Frtered) 


REPORT  DOCUMENTATION  PAGE ME^cgSJSfS^. 


2.  GOVT  ACCESSION  NO.j  3.  «EC|P| 


:por t numbzr  _ . A / 

AFOSR  - TR-  7 5 " 0 2 0 6 


4.  title  (end  Subtm 9) 

application  of  advanced  methods  for  identificatio: 

AND  DETECTION  OF  NUCLEAR  EXPLOSIONS  FROM  THE 
ASIAN  CONTINENT 


S.  TYPE  OF  REPORT  & PERIOD  COVERS: 

Interim 


6.  performing  org.  report  number 


7.  AUTHOR!'*} 

T._  C.  Bache 
J._  T.  Cherry 
J.  M.  Savino 

8.  CONTRACT  OR  GRANT  NUMBER'S} 

F44620-74-C-0063 

s.  performing  organization  name  and  aooress 
Systems,  Science  & Software 
La  Jolla,  CA  94550 

10.  PROGRAM  ELEMENT,  PROJECT  TASK 
AREA  & WORK  UNIT  NUMBERS 

62701E 
A0  1827 

H.  CONTROLLING  OFFICE  NAME  AllO  AOORESS 

Advanced  Research  Projects  Agency 
1400  Wilson  Boulevard 
Arlington,  VA  22209 

12.  REPORT  OATE 

October  1974 

13.  number  of  pages 

MONITORING  AGENCY  NAME  A AOORESSf//  different  from  Controlling  Office) 

Air  Force  Office  of  Scientific  Research //(  p 
1400  Wilson  Boulevard 

IS.  SECURITY  CLASS,  (ot  this  report) 

Unclassified 

Arlington,  VA  22209 

1S».  OECLASSI  FI  CATION  LDWNGRADING 
SCHEOULE 

16.  DISTRIBUTION  STATEMENT  (ol  this  Report) 

Approved  for  public  release;  di«prihtiM«n  nr 

>1  -tm}  TpH  . 

1 

1 

1 

17.  DISTRIBUTION  STATEMENT  (ot  the  ebstrset  entered  In  Block  20,  II  dilferent  from  Report) 

18.  SUPPLEMENTARY  NOTES 

y 

19.  KEY  WOROS  (Continue  on  reverse  tide  II  necessmry  end  Identity  by  block  number) 

Reproduced  by 

NATIONAL  TECHNICAL 
INFORMATION  SERVICE 

U S Department  of  Commerce 
Springfield  VA  22151 

• 

20.  ABSTRACT  (Continue  on  reverse  side  if  necessary  end  Identity  by  block  number) 

This  report  presents  the  results  of  a continuing  effort  to  obtain  an  optimum 
multi-discriminant/detection  procedure  for  earthquakes  and  explosions  occurring 
within  the  Asian  Continent.  Significant  improvement  of  theoretical  predictions 
of  teleseismic  ground  motion  from  explosions  was  realized  using  a deterministic 

computer  model  which  merges  nonlinear  shock  code  calculations  ,,’ith  far-Iie1  ’ 
stress  v:ave  propagation  through  a realistic  earth  structure.  The  addition  to 

tl"  aismic  codes  of  the  Thomson-HaskeJ.l  matrix  formulation  and  a split  eartn 


UNCLAS 


sifiWES  subject  to  c 


SEC  11  wy~Y  CLASSIFICATION  of  THIS  PAOEfThwi  Data  Entcrtd) 


structure  capability  provides  for  the  inclusion  of  more  realistic  earth 
structure  models  in  the  near  source  and  receiver  regions  and  for  the  intervening 
path  in  the  upper  mantle.  The  resulting  theoretical  seismograms  show  a striking 
similartiy  to  observations  and  indicate  that  a successful  matched  filter  for 
event  detection  and  discrimination,  based  on  a substantial  portion  of  an  event 
signal,  can  be  constructed. 


)C^ 


J 


SECURITY  CLASSIFICATION  OF  THIS  PAGE'W/ion  Dnta  hntered) 


mm 


SYSTEMS.  SCIENCE  AND  SOFTWARE 


SSS-R-7  5-2483 


APPLICATION  OF  ADVANCED  METHODS  FOR  IDENTIFICATION  AND 
DETECTION  OF  NUCLEAR  EXPLOSIONS  FROM  THE  ASIAN  CONTINENT 


T.  C.  Bache 
J.  T.  Cherry 
J.  M.  Savino 


D D C 


Semi-Annual  Report 


npoPJlDjaEP 

FFB  7 1975 

liETSIEITiTE 


D 


Contract  ho.  F44620-74-C-0063 


S3  Project  261 


October  1 97 4c  ; rc  'r!  u.  ' 


m rcr.fr  orr.cf.  of  r«c:f-'’T,r,n  research  (afsc) 

KOTICE  0?  TRAui'I.VTTM.  TO  DOC 
V tv'  ri\  ff  art  It,  been  rev :v.*ed  and  is 

: r ;to  !AW  AFR  130-12  (7b). 


[ •'  i d. 

D.  VL  TAY10R 

Tec'iii.cnl  Ir.Sorm-tTon  Onicer 


P O BOX  16. !0,  LA  JOLLA.  CALIFORNIA  92037.  TELEPHONE  (714)  453  0060 


. a 


R-2483 


FOREWORD 


This  semiannual  technical  report  entitled,  "Application 
of  Advanced  Methods  for  Identification  and  Detection  of  Nuclear 
Explosions  from  the  Asian  Continent,"  is  submitted  by  Systems, 
Science  and  Software  (S3)  to  the  Advanced  Research  Projects 
Agency  and  to  the  Air  Force  Office  of  Scientific  Research 
(AFOSR) . This  report  presents  the  results  of  a continuing 
effort  to  obtain  an  optimum  multi -discriminant/detection  pro- 
cedure for  earthquakes  and  explosions  occurring  within  the 
Asian  Continent.  The  work  is  being  performed  under  Contract 
Number  F44620-74-C-0063.  Mr.  William  J.  Best  is  the  AFOSR 
technical  contracting  officer. 

Dr.  J.  Theodore  Cherry  is  the  S3  project  manager.  Dr. 
Thomas  C.  Bache  is  responsible  for  the  development  and  appli- 
cation of  the  seismic  ground  motion  prediction  work.  Dr.  John 
M.  Savino  is  responsible  for  the  analysis  of  the  seismic  data 
against  which  all  theoretical  developments  must  eventually  be 
tested.  Acting  as  consultants  on  the  project  are  Professors 
Charles  B.  Archambeau  of  the  University  of  Colorado  and  David 
C.  Harkrider  of  the  California  Institute  of  Technology. 


ABSTRACT 


Significant  improvement  of  theoretical  predictions  of 
teleseismic  ground  motion  from  explosions  was  realized  using 
a deterministic  computer  model  which  merges  nonlinear  shock 
code  calculations  with  far-field  stress  wave  propagation 
through  a realistic  earth  structure.  The  addition  to  the  seis- 
mic codes  of  the  Thomson-Haskell  matrix  formulation  and  a 
split  earth  structure  capability  provides  for  the  inclusion 
of  more  realistic  earth  structure  models  in  the  near  source 
and  receiver  regions  and  for  the  intervening  path  in  the  up- 
per mantle.  The  resulting  theoretical  seismograms  show  a 
striking  similarity  to  observations  and  indicate  that  a suc- 
cessful matched  filter  for  event  detection  and  discrimination, 
based  on  a substantial  portion  of  an  event  signal,  can  be  con- 
structed. 

Parameter  studies  for  different  near  source  and  re- 
ceiver earth  structures  clearly  demonstrate  how  the  conventional 
procedure  employed  for  body-wave  magnitude  calculations  intro- 
duces scatter  into  the  determinations.  A more  consistent 
technique  for  mfe  determinations  is  presented  here  and  is  ex- 
pected to  result  in  more  accurate  scaling  relationships  for 
events  of  different  yields. 

A variable  frequency  magnitude  technique,  designed  to 
exploit  spectral  differences  between  earthquakes  and  explo- 
sions, was  developed  and  applied  to  short-period  recordings 
of  a large  population  of  Eurasian  events.  Complete  discimina- 
tion  between  all  explosions  and  shallow  earthquakes  down  to 
m (LASA)  ^4.5  was  achieved.  This  technique  was  found  to 
discriminate  various  multiple  explosion  scenarios  as  well. 
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I.  INTRODUCTION 


O 


c 


0 


t; 


The  objective  of  the  research  program  described  in  ’is 
report  is  to  merge  relevant  knowledge  from  rock  mechanics, 
numerical  simulation  of  stress  wave  propagation,  information 
theory  and  seismology  to  obtain  an  optimum  multi-discriminant/ 
detection  procedure  for  earthquakes  and  explosions  occurring 
within  the  Asian  continent.  During  the  first  hal  year  o e 
program,  work  has  proceeded  in  several  task  areas  including: 

1.  Refine  and  improve  the  S3  capability  for  a deter 
ministic  prediction  of  teleseismic  ground  motion 
due  to  underground  explosions  and  earthquakes; 

2.  Collect  and  analyze  ground  motion  data  from  a 
variety  of  seismic  networks; 

3.  Design  a multi-discriminant  procedure  for  detect- 
ing and  discriminating  nuclear  explosion  events; 

4.  Test  suggested  detection/discrimination  procedures, 
both  conventional  and  newly  designed,  using  actual 
seismic  data. 

The  seismic  code  package  now  operational  at  S1  consists 
of  more  then  twenty  separate  computer  codes  written  over  a 
period  of  years  by  many  different  programmers.  These  codes  in 
elude  the  capability  to  compute  the  elastic  wave  generate  y 
a realistic,  large  displacement,  nonlinear  explosion  source  an 
to  propagate  the  resulting  body  and  surface  waves  to  teleseis 
mic  distances  where  synthetic  seismograms  are  gen  ated. 
code  package  also  includes  a versatile  seismic  data  ana  ysis 

capability. 

Recent  additions  and  improvements  to  the  seismic  codes 
are  discussed  in  Section  II  and  III  of  this  repoit.  Of  par- 
ticular interest  is  the  addition  of  crustal  reverberations 
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to  the  body  wave  propagation  codes.  Preliminary  results  illus- 
trating the  influence  of  crustal  structure  on  the  seismograph 
record  are  presented  in  Section  II. 

A recently  suggested  variable  frequency  magnitude  dis- 
criminant^22^ employing  spectral  magnitudes  computed  from  short 
period  compressional  waves  is  being  carefully  tested  and  is 
discussed  in  Section  IV. 

Finally,  the  application  of  our  modeling  capability 
to  study  the  dependence  of  seismic  magnitudes  on  yield,  near 
source  material  properties,  depth  of  burial  and  earth  struc- 
ture was  the  subject  of  an  S3  report  prepared  in  August  at 
AFOSR  request This  report  is  included  as  Appendix  A. 
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II.  DETERMINISTIC  PREDICTION  OF  TELESEISM  I C GROUND  MOTION 
FROM  EARTHQUAKES  AND  UNDERGROUND  EXPLOSIONS 
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Combining  relevant  knowledge  from  rock  mechanics, 
numerical  simulation  of  stress  wave  propagation  and  seismo- 
logy, a deterministic  prediction  of  teleseismic  ground  motion 
from  earthquakes  and  underground  nuclear  explosions  has  been 
developed  at  S3  and  is  rapidly  being  improved.  This  capabi- 
lity and  its  application  to  a variety  of  problems  is  des- 
cribed in  References  1-4  and  7. 

Fundamental  to  the  computational  techniques  applied 
is  the  representation  of  the  source  in  terms  of  an  equivalent 
elastic  source.  Various  linear  wave  propagation  techniques 
may  then  be  applied  to  compute  teleseismic  body  and  surface 
waves.  The  equivalent  elastic  source  may  be  obtained  for  a 
variety  of  earthquake  and  explosion  source  models  including: 

1.  A realistic  large  displacement,  nonlinear  finite 
difference  calculation  of  the  explosion  generated 
shock  wave  and  its  conversion  to  an  elastic  wave. 
In  most  cases  a one -dimensional  (spherically 
symmetric)  calculation  will  suffice;  when  re- 
quired, two-dimensional  (axisymmetric)  calcula- 
tions have  also  been  used. 

2.  The  relaxation  of  the  prestressed  material  sur- 
rounding the  explosion  induced  cavity  (tectonic 
release)  using  a theoretical  model  developed  by 
Archambeau . ^ 

3.  An  analyti.:al/numerical  earthquake  model  of 
the  relaxation  type  developed  by  Archambeau 
and  Minster.  ^ 

Having  the  equivalent  elastic  source  representation, 
Rayleigh  and  Love  surface  waves  are  computed  using  a 
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computational  technique  developed  by  Harkrider1  1 based  on 
the  Thomson -Haskell  matrix  method. 

The  propagation  of  body  waves  is  accomplished  using 
a combination  of  geometric  ray  theory  and  Thomson -Haskell 
matrix  methods  in  a layered  earth  model.  The  anelastic 
attenuation  parameter,  generally  denoted  by  the  material 
quality  factor,  Q,  is  included  in  this  analysis. 

The  remainder  of  this  section  includes  a more  de- 
tailed description  of  the  equivalent  elastic  source  repre- 
sentation and  the  techniques  used  to  compute  teleseismic 
body  and  surface  waves,  in  particular  those  features  to 

r 4i 

which  improvements  have  been  made  since  our  August  report1  1 
which  is  included  as  Appendix  A. 


2.1  EQUIVALENT  ELASTIC  SOURCE  REPRESENTATION 


The  Fourier  transformed  equations  of  motion  in  a 
homogeneous,  isotropic,  linearly  elastic  medium  may  be 
written 


where  u is  particle  displacement  and  k and  k are  the 

p s 

compress ional  and  shear  wave  numbers.  The  Cartesian  poten- 
tials x and  X are  defined  by 

4 


(2.2) 


and  may  be  easily  shown  to  satisfy  the  wave  equation 


0,  a = 1,2, 3, 4, 


(2.3) 
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o)/vs,  i = 1,2,3. 


This 


where  k = k = w/v  , k-  = k 
4 P piS 

equation  has  as  a solution  the  following  expansion  in 
spherical  eigenfunctions  (e.g.,  Morse  and  Feshbach , ) , 


Xa(r,uO 


X)  hi2)tkor)  S t")  C0Sm*  * 5inm<1' 

n=0  m*0 


*Pm(cos0) 

n 


(2.4) 


where  the  h^^  are  spherical  Hankel  functions  of  the  second 


n 


kind  and  the  Pm  are  associated  Legendre  functions.  The 
n 

vector  r has  as  components  the  spherical  coordinates  r,0,4>. 

Equations  (2.4)  together  with  (2.1)  provide  an  equiva- 
lent elastic  source  representation  of  the  outgoing  elastic 
waves.  The  values  of  the  multipole  coefficients  (u)  , 

(u)  , a = 1,2, 3, 4,  prescribe  the  displacement  field  at  all 
points  in  the  homogeneous  elastic  medium  where  (2.1)  applies. 

Given  the  displacement,  u^,  or,  alternatively,  the 
potentials  xa>  °ne  may  determine  the  multipole  coefficients 
using  (2.4).  For  example,  using  the  orthogonality  of  the 
spherical  eigenfunctions,  one  may  derive 


A(“>M 

nm 

nm  k ' 
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That  is,  given  the  displacements,  u(r,w)  on  a sphere  of 
radius  r,  one  can  use  (2.2),  (2.5)  and  (2.4)  to  compute  the 
displacements  at  any  point  in  the  elastic  regime.  This  pro- 
cedure was  first  suggested  for  geophysical  applications  by 
Archambeau. ^ ^ 

Using  the  equivalent  elastic  source  representation, 
highly  nonlinear  finite  difference  calculations  of  the  near 
source  earthquake  or  explosion  displacement  fields  may  be 
coupled  to  the  elastic  field  where  efficient  elastic  wave 
propagation  techniques  may  be  applied.  Carrying  the  nonlinear 
computation  out  to  a radius  at  which  the  material  response 
is  in  the  elastic  regime,  Eq.  (2.5)  may  be  applied  to  compute 
the  multipole  coefficients  and  thereby  the  equivalent  elastic 
source . 

For  a spherically  symmetric  explosion,  the  explosion 
induced  elastic  field  is  often  represented  by  the  reduced 
displacement  potential,  t^(x),  defined  by  (e.g.,  Cherry, 
et  al . ) , ^ 


ur(r,t) 


_i_  Mil 
rv  dr  ’ 
P 


(2.6) 


where  x * t - r/Vp,  is  the  retarded  time.  Note  that  <Kx) 
is  independent  of  x. 


Applying  the  Fourier  transform,  (2.6)  becomes 

(1  + ik  r)  -ik  r „ 
u (r,w)  * 2 — e p 4>(u>)  . 

1 -r-2 


(2.7) 


Then,  carrying  out  the  algebra  (observing  that  the  u of 
(2.2)  is  the  Cartesian  displacement  vector),  it  may  be  shown 
that  the  multipole  coefficients  reduce  to  a single  nonzero 
term,  the  monopole: 
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a(4) 


0 0 


(«) 


• ikp  <Kw) 


(2.8) 
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For  more  complicated  sources,  additional  nonzero  co- 
efficients are  present.  For  example,  in  the  case  of  an  ex- 

^(4) 

plosion  in  an  axisymmetric  cavity  the  coefficients  , 

g(l),  a^2^,  n = 2,4,6  ...  , are  also  nonzero. ^ Also, 

n » n i r a i 

Archambeau1  J has  given  a theory  for  tectonic  stress  release 
due  to  underground  explosions  which  gives  nonzero  quadrupole 


coefficients  A 


a(«) 


2m 


B 


(°0 

2m 


m 


0,1,2. 


2.2  TELESEISMIC  BODY  WAVES 


The  basic  technique  used  for  propagating  source  (given 
by  an  equivalent  elastic  source  representation)  generated 
body  waves  through  a layered  earth  structure  is  geometric 
ray  theory  in  a layered  earth  model,  including  the  anelastic 
attenuation  parameter  Q.^  One  layered  earth  model  may 
be  used  to  represent  the  earth  structure  beneath  the  source 
while  a lecond  represents  the  earth  structure  in  the  receiver 
half  of  the  travel  path.  The  ray  crosses  into  the  second 
(receiver)  structure  at  its  point  of  deepest  penetration 
which,  for  teleseismic  distances,  is  deep  within  the  mantle 
where  lateral  inhomogeneities  are  thought  to  be  small. 

The  unembellished  ray  theory  described  above  has  been 

• quite  successful  in  illuminating  the  important  gross  proper- 
ties of  the  teleseismic  body  wave  record  (see  Appendix  A). 
Additional  complexities  observed  on  the  records  are  thought 
to  be  the  result  of  crustal  structure  details  in  the  source 

t and  receiver  regions.  These  details  may  be  included  using 

Thomson-Haskell  matrix  techniques. 

The  propagation  of  body  waves  is  shown  schematically 
in  Fig.  2.1.  Referring  to  the  figure,  the  options  available 

* for  propagating  body  waves  are  as  follows: 
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Figure  2.1.  Schematic  view  of  body  wave  propagation  from  a source  at  a depth 
a receiver  at  teieseismio  distances.  Four  layered  earth  models  may  be  used  ti 
sent  different  portions  of  the  path:  (1)  Crustal  structure  of  thickness  Ds 

source  region;  (2)  Mantle  structure  for  the  source  region;  (3)  Mantle  structu 
the  receiver  region;  (4)  Crustal  structure  of  thickness  Dn  for  the  receiver 
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1.  Geometric  ray  theory  with  an  average  earth 
structure  for  the  entire  travel  path 
(available  in  1973) . 

2.  Geometric  ray  theory  with  two  earth  structure 
models,  one  for  the  source  region  and  one  for 
the  receiver  region  (used  in  preparing  Appendix 
A in  August)  . 

3.  Use  Haskell  matrices  for  receiver  crustal 
structure  (recently  added). 

4.  Use  Haskell  matrices  for  source  region  crustal 
structure  (recently  added) . 

Therefore,  if  us(w)  is  the  far-field  component  of 
the  source  generated  displacement  spectrum,  then  the  ground 
motion  spectrum,  u^(w) , at  teleseismic  distances  is 

uR(«)  - us(w)  • 1^(u)  • t£(u))  • Tg  (w)  • T“(w)  , (2.9) 

where  TG  and  TG  represent  transfer  functions  computed  by 
S K 

ray  theory  for  the  source  and  receiver  structures  and 

T**  ?nd  Tn  represent  transfer  functions  computed  by 
s ** 

Haskell  matrices  for  the  source  and  receiver  structures. 

The  ray  theory  transfer  functions  include  attenuation  and 
dispersion  via  the  material  quality  factor,  Q.  The  effec- 
tive Q for  the  ray  travel  path  is  used  in  the  multiplica- 
tive factor 


where  T is  the  travel  time  and  is  a high  frequency 

cutoff,  generally  set  to  * 2000n  (the  results  are  rather 
insensitive  to  this  parameter). 
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The  appropriate  theory  for  the  transfer  function  TnO) 

r I t 1 K 

has  been  given  by  Haskell.1  J Using  his  notation,  for 
incident  P '•'aves  the  transfer  function  is 


tr(u)  =■  2c(J^  - Jsi)/an  D , (2.11) 

for  the  vertical  component  and 

Tr(«)  = 2cCj,2  • J32^/an  D ’ C2  *12) 


for  the  radial  component  of  motion.  The  phase  velocity 
c * an  c sc C 4> where  is  computed  by  ray  theory  (Fig. 

2.1). 

LI 

At  present  the  transfer  function  Tg(u>)  can  only  be 
computed  for  the  case  of  the  spherically  symmetric  point 
source  for  which  the  equivalent  elastic  source  is  given  by 
Eq.  (2.8).  The  transfer  function  for  the  wave  emerging  at 
an  angle  $s  from  the  source  crustal  structure  (see  Fig.  2.1), 
including  a point  source  at  a depth  of  burial  H,  has  been 
derived  by  Fuchs. Normalizing  on  the  displacement  for 
the  same  source  in  a halfspace  of  the  source  layer  material, 
this  transfer  function  may  be  written 


TsM 


2 ps  sin  4>s  cos 


+ 


(2.13) 


where  the  notation  of  Fuchs,  Eq.  (3.26),  has  been  used. 

Fuchs  also  gives  the  appropriate  transfer  function 
for  a double  couple  source  which  can  be  modified  to  compute 
T^(w)  for  an  explosion  plus  tectonic  release.  However,  this 
has  not  yet  been  implemented. 

It  would  be  desirable  to  have  a general  theory  by 
which  an  equivalent  elastic  source  including  any  number  of 


10 


ffH». I HWIW  "I I' H»JWW 


R-  2483 


* 


* 


I 


multipole  coefficients  could  be  introduced  into  a Haskell 
matrix  representation  of  the  crustal  structure.  Such  a 
theory  has  been  developed  for  surface  waves  by  Harknder 
and  Archambeau,[1S]  setting  forth  the  procedure  to  be 

followed. 

At  present,  the  attenuation  parameter  Q is  not 
included  in  the  Haskell  matrix  computations.  However,  an 
elementary  modification  by  which  inclusion  of  Q may  be 
accomplished  has  been  suggested  by  Hasegawa.  This 

will  be  implemented  in  the  near  future. 


2.3  TE LESEISMIC  SURFACE  WAVES 


Consider  a source  of  pure  dilatation  m one  layer, 

denoted  s,  of  a multilayered  material.  Then^the  spectrum 

of  the  vertical  Kayleigh  wave  displacement,  w , is  given 
r q i 

by  Harkrider1  J as 


-i4irysvs(  )kRKR*R  e 


IT. 


Ht2’(kRr) 


(2.14) 


■ 

I 

r. 


i 


where 


y * rigidity  of  the  material  in  layer  s, 
ij)  (w)  * Fourier  transformed  reduced  displacement 
potential , 

k^  = oj/cj^,  where  is  the  Rayleigh  wave 

velocity  at  the  frequency  uy 

K = depth  dependent  Rayleigh  wave  excitation 
K 

factor , 

A = amplitude  response  of  the  multilayered 
— R 

material , 

t - uR  /a  , where  R is  the  effective 
s o s 0 . 

cavity  radius  for  the  explosion  and 
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as  is  the  P wave  velocity  in  layer 
s, 

H(2)(kDr)  ■ spherical  Hankel  function  of  the 
0 K 

second  kind. 


A similar  formula  applies  for  a source  specified  by  a general 
equivalent  elastic  source  representation. 

Surface  wave  calculations  using  this  formulation  are 
presented  in,  for  example,  Appendix  A.  However,  one  weakness 
of  this  approach  is  that  only  a single  (average)  earth 
structure  model  can  be  used  for  the  entire  source-receiver 
travel  path.  For  surface  waves  traversing  significant 
discontinuities  such  as  those  at  oceanic-continental  bound- 
aries, large  errors  may  be  introduced. 

In  order  to  improve  our  surface  wave  computational 
capability,  our  computer  programs  are  presently  - 'ing 
modified  to  include  the  option  of  propagation  over  a mixed 
path.  In  particular,  the  path  may  be  composed  of  two  earth 
structure  models;  the  source  medium,  1,  and  the  receiver 
medium,  2.  While  the  modifications  will  be  valid  for  a 
general  source,  it  is  once  again  convenient  to  present  them 
in  the  context  of  the  center  of  dilatation. 


For  a mixed  path,  the  approximate  vertical  displace- 
ment Rayleigh  wave  spectrum  for  a selected  mode  is  given  by 

rH,  ! 

T IV  ^ 1 | c 

X 1 


{w  } 
0 


•^Sl^S,^ 


kRiKR,^RxTRi2Ho 


(2)  (kR  r .+K 


R 2 r 2 ^ 6 


(YD.r  +yR  r ) 


L2  2 


(2.15) 


where  the  subscript  denotes  media  1 or  2.  The  total  path 
r is  defined  as 


r ■ r + r 
1 2 


(2.16) 
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Assuming  that  the  total  horizontal  energy  flux  for  a 
given  Rayleigh  wave  mode  remains  constant  in  the  process  oi 
transmission  across  the  boundary  separating  the  two  media, 
we  have 


(2.17) 


where  WRi  is  the  total  horizontal  energy  flux  in  medium  i 
normalized  to  a unit  vertical  displacement  at  the  surface. 

r 171 

The  normalized  energy  flux  is  given  by1  J 


O 


(2.18) 


where  is  the  group  velocity  and  is  the  total  potential 

or  kinetic  energy  at  a particular  frequency  normalized  to  the 
surface  vertical  displacement. 

This  correction  factor,  which  has  given  satisfactory 
results  in  predicting  the  amplitude  changes  of  plane  Rayleigh 
waves  traveling  across  a lateral  boundary,  has  been  given  by 
Q McGarr^18^  and  McGarr  and  AlsopJ19^  It  assumes  that  there 

is  little  or  no  mode  conversion  or  reflection  at  the  boundary, 
that  the  boundary  is  far  enough  from  the  source  so  that  the 
wave  front  has  small  lateral  curvature  and  that  the  wave  front 
;;  is  parallel  to  the  boundary,  i.e.,  no  refraction  effects  are 

included. 

The  validity  of  this  approximation  has  been  verified 
by  wave  propagation  studies  in  laterally  heterogeneous  media 
f using  the  finite  difference  numerical  technique  (D.  Boore, 

personal  communication  to  R.  Alewine)  . 

For  more  than  two  media  in  the  propagation  path  we  can 
generalize  under  the  same  assumption  to 
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1 


r 

R‘n  \EnUn/ 


(2.19) 


and  the  propagation  phase  and  decay  terms  become 


h'2)(  E Wi) 

\ i-i  ' 


(2.20) 


E Vi 

i-1 


(2.21) 


respectively. 

The  implementation  of  the  theory  described  above  is 
expected  to  substantially  enhance  our  surface  wave  computa- 
tional capability. 
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2.4  EFFECTS  OF  CRUSTAL  STRUCTURE  ON  EXPLOSION  GENERATED 

V 'WAVES' 

Preliminary  results  illustrating  the  effect  of  the 
crustal  structure  in  the  source  and  receiver  regions  are 
presented  in  this  section. 

Consider  a typical  explosion  of  moderate  yield.  The 
far-field  component  of  the  displacement  pulse  generated  by 
such  an  explosion  is  shown  in  Fig.  2.2.  Plotted  is  <Kt)/v 
which,  at  large  r,  is  the  far-field  displacement  (see 
Eq.  (2.6)).  The  time  scale  is  appropriate  to  a yield  of 
20  kT. 

According  to  Eq.  (2.9),  propagation  of  the  source 
generated  stress  wave  to  teleseismic  distances  is  accom- 
plished by  convoluting  with  transfer  functions  representing 
the  effects  of  the  crust  and  upper  mantle.  Synthetic  seismo- 
grams can  then  be  generated  by  convolution  with  a transfer 
function  representing  the  response  of  the  seismometer  of 
interest . 

Consider  first  the  synthetic  seismogram  computed  using 
geometric  ray  theory  to  represent  the  entire  travel  path. 

In  Fig.  2.3  is  shown  such  a seismogram  for  an  explosion  at  a 
burial  depth  such  that  the  average  source- to-surface  travel 
time  for  P waves  is  0.2  sec.  The  epicentral  distance  is 
approximately  4000  km  and  the  ray  path  is  characterized  by 
T/Q  = 1.1,  where  T is  the  travel  time.  Since  the  arrival 
time  lag  between  the  direct  (P)  and  free-surface  reflected 
(pP)  phases  is  so  small  (=  0.4  sec),  the  wave  train  appear- 
ing on  the  seismogram  represents  a combination  of  the  two 
phases.  The  amplitudes  given  on  this  and  succeeding  seismo- 
grams are  in  millimicrons.  Also,  the  time  scale  on  this 
and  succeeding  seismograms  is  relative  in  that  the  zero  point 
has  no  consistent  relation  to  explosion  detonation  time. 
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The  nature  of  the  effect  of  the  source  region  crustal 

* structure  on  the  seismogram  may  be  seen  by  examining  a series 
of  seismograms  for  which  the  transfer  function  due  to  the 
source  crustal  structure  has  been  included  in  the  calculation. 
The  P wave  velocity  profiles  for  several  crustal  structures 

^ typical  of  NTS  are  shown  in  Fig.  2.4.  As  seen  in  the  figure, 

the  main  differences  between  these  structures  are  restricted 
to  the  top  kilometer.  In  Fig.  2.5  synthetic  seismograms 
illustrating  the  effect  of  the  source  crustal  structure  SI 

* of  Fig.  2.4  are  shown.  In  Fig.  2.5a  the  upper  mantle  model 
used  for  the  ray  theory  portion  of  the  calculation  was 
identical  to  that  used  for  computing  the  seismogram  of 
Fig.  2.3.  Therefore,  comparing  Figs.  2.3  and  2.5a,  the  ef- 

* feet  of  adding  to  the  calculation  a detailed  model  of  the 
source  region  crustal  structure  is  seen.  When  making  this 
comparison,  it  should  be  noted  that  two  different  instrument 
response  functions  were  used  for  the  two  calculations.  How- 
ever, as  will  be  demonstrated  below,  this  has  a minor  effect. 

Two  additional  seismograms,  computed  using  the  source 
crustal  structure  SI,  are  shown  in  Figs.  2.5b  and  2.5c.  For 
these  two  seismograms  a higher  Q was  introduced  into  the 
upper  mantle  model  for  the  ray  theory  portion  of  the  calcula- 
tion resulting  in  T/Q  = 0.7.  The  effect  of  the  two  instrument 
response  functions,  denoted  #1  and  #2,  are  seen  by  comparing 
Figs.  2.5b  and  2.5c.  Finally,  in  Fig.  2.5d  the  effect  of  an 
intermediate  Q model,  T/Q  a 0.84,  is  demonstrated. 

The  seismograms  of  Fig.  2.5  were  computed  using  a 
single  model  of  the  source  crustal  structure.  This  model, 

/ denoted  SI  in  Fig.  2.4,  is  characterized  by  a low  velocity 

zone  with  a high  velocity  cap  in  the  region  between  source 
and  free  surface.  A few  hundred  meters  below  the  source  is 
a sharp  velocity  contrast  at  the  Paleozoic  rock  which  under- 
lies the  tuff.  The  remainder  of  the  crustal  structure  down 
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to  the  M discontinuity  is  taken  from  Hill  and  Pakiser. 

The  effect  of  this  crustal  structure  is  to  give  the  seismo- 
gram its  commonly  observed  "W"  shape  and  to  spread  its 
characteristic  period. 

The  effect  of  variations  in  the  uppermost  part  of  the 
source  crustal  structure  may  be  demonstrated  by  several 
examples.  For  the  seismogram  of  Fig.  2.6a  the  crustal 
structure  S2  was  used,  thus  bringing  the  Paleozoic  layer 
nearer  the  source.  Comparing  Figs.  2.6a  and  2.5a,  one  ob- 
serves that  the  Paleozoic  depth  primarily  influences  the 
amplitude  of  the  first  peak.  A similar  comparison  is  af- 
forded by  Figs.  2.6b  and  2.5b,  this  time  for  a different 
instrument,  with  the  main  effect  being  a change  in  the 
relative  amplitudes  of  the  first  and  second  troughs. 

The  effect  of  the  layering  above  the  source  is  por- 
trayed in  Fig.  2.7.  The  crustal  model  S4  gives  an  average 
velocity  profile  for  the  region  above  the  source,  smoothing 
the  high  and  low  velocity  detail  cf  models  SI,  S2,  and  S3. 

The  model  S5  uses  average  properties  all  the  way  down  to  the 
Paleozoic  rock.  Comparing  Figs.  2.7a  and  2.7b  to  2.5a,  it 
is  cltar  that  the  first  few  swings  of  the  seismogram  are 
strongly  influenced  by  the  detailed  structure  in  the  source 
region. 

The  crustal  model  S6  is  similar  to  S5  except  that  the 
M discontinuity  is  at  a much  shallower  depth.  The  resulting 
seismogram  is  shown  in  Fig.  2.8a  and  should  be  compared  to 
Figs.  2.5a  and  2.7b.  Finally,  crustal  model  S7  has  the  M 
discontinuity  essentially  deleted.  The  effect  on  the  pulse 
shape  is  shown  in  Fig.  2.8b. 

It  has  been  demonstrated  that  the  near-source  crustal 
structure  can  have  a strong  influence  on  the  seismograph 
trace.  The  effect  of  a few  typical  receiver  regional  crustal 
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a)  Source  crustal  model  S2,  Instrument  #2,  T/Q  * 1.1. 


*IM  UfCl 

b)  Source  crustal  model  S3,  Instrument  #1,  T/Q  * 0.7. 

Figure  2.6.  The  effect  of  a modified  source  crustal  structure 
with  the  underlying  Paleozoic  rock  at  shallow  depth. 
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structures  on  the  wave  form  will  now  be  demonstrated. 

Velocity  profiles  for  the  receiver  region  crustal  structures 
R1  and  R2  are  shown  in  Fig.  2.9.  The  only  difference  between 
the  two  is  the  properties  of  the  top  layer  — one  structure 
corresponds  to  an  instrument  on  bedrock  and  the  other  to  an 
instrument  in  the  same  area  located  on  the  sediments. 

In  Fig.  2.10a  is  shown  the  seismogram  including  the 
effects  of  both  source  structure  SI  and  receiver  structure 
Rl.  This  should  be  compared  to  Fig.  2.5a  where  only  SI  was 
included.  A new  parameter  is  varied  for  Fig.  2.10b  which 
was  calculated  with  a different  depth  of  burial.  In  Fig. 
2.10c  the  structure  models  S3  and  Rl  are  used  (compare  to 
Fig.  2.6b).  Finally,  in  Fig.  2.11  three  further  examples 
of  various  combinations  of  source  crustal  structure  models 
and  the  receiver  model  R2  are  shown. 

In  Table  2.1  the  important  parameters  used  in  the 
calculations  leading  to  the  seismograms  presented  in  this 
section  are  summarized  along  with  key  measurements  taken 
from  the  seismograms.  The  quantities  A^,  A^,  A^,  T,  are 
defined  in  Fig.  2.11c. 

Finally,  a brief  description  of  the  distinguishing 
features  of  the  various  crustal  models  of  this  section 
(Figs.  2.4  and  2.9)  is  given  in  Table  2.2. 
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c)  DOB  ■ 0.4  km,  Source  crustai""structure  S3,  T/Q  = 0.7. 


Figure  2.10.  Synthetic  seismograms  including  both  the  in- 
dicated source  crustal  structures  and  the  receiver  crustal 
structure  Rl.  Instrument  If 2 was  used. 
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a)  Source  crustal  structure  SI,  Instrument  #2,  T/Q 
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b)  Source  crustal  structure  S8,  Instrument  #1,  T/Q 
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c)  Source  crustal  structure  S3,  Instrument  #1,  T/Q  ■ 0.84. 

Figure  2.11.  Synthetic  seismograms  illustrating  the  effect 
of  several  source  crustal  structure  models  in  combination 
with  receiver  crustal  structure  model  R2. 
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TABLE  2.1 


SUMMARY  OF  THE  IMPORTANT  PARAMETERS  CHARACTERIZING 
THE  SEISMOGRAMS  OF  FIGURES  2.3  - 2.11 


* 

Source 

Receiver 

A 

A 

A 

Figure 

DOB 

Structure 

Structure 

Instr. 

T/Q 

_X- 

T 

2.3 

0.4 

None 

None 

#1 

1.1 

7.7 

5.7 

0.55 

1.30 

2.5a 

0.4 

SI 

None 

»2 

1.1 

5.7 

7.1 

5.2 

1.10 

z 

2.5b 

0.4 

SI 

None 

#1 

0.7 

6.3 

8.8 

6.1 

, . J2 

2.5c 

0.4 

SI 

None 

#2 

0.7 

5.8 

8.6 

6.9 

1.05 

2 . 5d 

0.4 

SI 

None 

#2 

0.84 

6.0 

8.2 

6.1 

1.05 

T» 

%» 

2.6a 

0.4 

S2 

None 

#2 

1.1 

5.7 

7.3 

5.4 

0.98 

2.6b 

0.4 

S3 

None 

#1 

0.7 

6.2 

8.8 

6.5 

0.95 

2.7a 

0.4 

S4 

None 

12 

1.1 

6.2 

6.8 

4.4 

1.01 

2.7b 

0 . 5* 

SS 

None 

*2 

1.1 

6.2 

5.6 

3.6 

1.32 

•* 

2.8a 

0 . 5* 

S6 

None 

*2 

' 1.1 

6.5 

7.4 

5.0 

1.15 

2.8b 

0.4 

S7 

None 

12 

0.7 

7.5 

10.7 

8.0 

1.05 

2.10a 

0.4 

SI 

R1 

#2 

1.1 

8.5 

10.7 

7.5 

1.05 

z 

2.10b 

0.5 

si 

R1 

*2 

1.1 

8.6 

11.7 

8.6 

1.10 

2.10c 

0.4 

S3 

R1 

*2 

0.7 

8.5 

12.7 

10.7 

0.97 

2.11a 

0.4 

SI 

R2 

*2 

0.7 

7.0 

10.1 

8.2 

1.05 

2.11b 

0.4 

S8 

R2 

*1 

0.7 

7.9 

10.4 

7.0 

0.97 

t 

2.11c 

0.4 

S3 

R2 

#1 

0.84 

7.5 

9.8 

7.3 

9.95 

A depth  of  0.5  km  is  required  in  source  crustal  structures  S5  and  S6  to 
maintain  the  source -surface  travel  time  at  a constant  value  for  all  seismo- 
grams of  the  table  excepting  2.10b. 
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Crustal 

51 

52 

53 

54 

55 

56 

57 

58 

R1 

R2 


TABLE  2.2 

SUMMARY  OF  THE  MAJOR  FEATURES  OF  THE  CRUSTAL 
MODELS  SI  - S8  and  R1  - R2 

Model  Important  Features 

Standard  source  crustal  model  with 
Paleozoic  rock  at  0.7  km  and  M 
discontinuity  at  31  km 

Paleozoic  rock  at  0.45  km 

Paleozoic  rock  at  0.475  km 

Average  properties  in  top  0.35  km 

Average  properties  in  top  0.98  km 

Average  properties  in  top  0.7  km, 

M discontinuity  at  18.5  km 

M discontinuity  not  present 

SI  modified  to  reduce  velocity 
contrast  in  top  0.5  km 

Standard  receiver  crustal  model  with 
sedimentary  top  layer 

R1  modified  for  bedrock  top  layer 


32 


R- 2485 


III.  DATA  ACQUISITION,  EDITING  AND  ANALYSIS 

3.1  DATA  ACQUISITION  AND  EDITING 

During  the  first  half-year  of  the  contract  period,  con- 
siderable time  was  spent  quiring  and  editing  numerous  seis- 
mic data  sets  for  North  American  and  Eurasian  earthquakes  and 
explosions.  These  data  form  the  basis  for  testing  our  wave- 
form modeling  capabilities,  the  applicability  of  newly  pro- 
posed discriminants  (e.g.,  see  Section  IV  on  Spectral  Magni- 
tude Discriminant)  and  the  specification  of  a matched  filter 
for  improved  event  detection  and  identification. 

At  the  present  time  our  library  of  Eurasian  events  con- 
sists of  192  earthquakes  and  42  presumed  explosions  recorded 
on  digital  magnetic  tape.  One  hundred  ninety  of  these  event 
signatures  are  short-period  recordings  of  best  beam  signals 
from  the  full  LASA  array  while  the  remaining  44  events  are 
best  beam  short-period  signals  from  the  original  Oyer  subarray 
at  Norway.  These  data  have  been  acquired  through  the  coopera- 
tion of  people  at  Lincoln  Laboratory,  in  particular,  Dr.  Richard 
Lacoss . 

Digital  magnetic  tape  data,  both  short-  and  long-period, 
from  the  former  LRSM  stations,  have  been  obtained  from  SDL  for 
a number  of  NTS  explosions.  These  data  represent  multi-station 
recordings  of  individual  explosions  and  provide  important  in- 
formation on  different  combinations  of  source-to-receiver  path 
conditions.  A request  submitted  to  SDL  for  additional  events 
is  pending. 

Each  of  the  sets  of  digital  data  has  required  the  de- 
velopment of  an  editing  program  for  unpacking  and  verifying 
the  data.  As  part  of  the  editing  procedure  seismograms  were 
plotted  on  a conventional  time  base  for  all  of  the  above  events 
and  comprise  a very  important  part  of  our  event  library. 
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In  addition  to  the  above  digitally  formatted  data, 
copies  of  seismograms  from  selected  high-gain  long-period 
stations  and  from  the  World-Wide  Network  of  Standard  Seismo- 
graph Stations  (WWNSS)  have  been  obtained.  These  data  pro- 
vide both  extended  geographical  coverage  and  a broader  band 
look  at  event  signals.  Our  plans  are  to  continue  expanding 
the  data  library  to  include  a wide  variety  of  event  and  re- 
ceiver settings,  so-called  anomalous  (Mg-m^  wise)  earthquakes, 
and  intraplate  earthquakes. 

3.2  DATA  ANALYSIS 

In  addition  to  the  extensive  code  capability  for  model- 
ing seismic  phenomena  developed  at  S3  and  described  in  Section 
II  of  this  report,  two  codes  were  recently  brought  on-line 
for  analysis  of  the  digitally  formatted  seismic  data.  The 
two  new  codes,  which  have  been  used  extensively  for  applica- 
tion of  the  spectral  magnitude  discriminant,  are  DATCON  and 
MAGIC.  DATCON  accepts  seismic  data  in  digital  format  and 
has  options  to  detrend,  taper,  demean  and  band-pass  filter 
the  data. 

MAGIC,  a MAGnitude  Identification  Code,  also  accepts 
seismograms  of  body  and/or  surface  waves  in  digitized  form. 

The  seismograms  may  be  real  data  or  synthetic  seismograms 
generated  by  source  codes  in  the  S3  library.  Several  events 
(explosions)  may  be  added  together  with  some  desired  spacing 
in  distance  and  time  to  give  a multiple  event  scenario.  The 
single  or  multiple  event  seismograms  may  then  be  plotted. 

A series  of  nariow  band  filters  with  center  frequen- 
cies, f,  and  quality  factors,  Q,  may  be  specified  for  recur- 
.ively  filtering  the  seismograms.  The  maximum  output  ampli- 
tude of  the  resulting  filtered  seismogram  is  then  picked  as 
the  amplitude  associated  with  the  energy  at  the  center  fre- 
quency of  the  filter.  Various  magnitude  definitions  may  then 
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be  applied  and  standard  Gutenberg-Richter  distance  corrections 
are  included  in  the  code.  The  code  generates  spectral  magni- 
tudes at  several  different  frequencies,  in^Cf),  over  time 
windows  of  varying  duration.  Plots  of  m^(f  ) versus  m^(f  ), 
f^  < f^,  can  then  be  tested  for  event  discrimination.  In 
addition  the  code  can  determine  and  fTg,  where  NTg  is  the 
spectral  surface  wave  magnitude  for  Rayleigh  waves  at  20 

second  period  and  is  the  same  magnitude  for  Love  waves, 

a — R 

and  test  for  event  discrimination  using  m,  versus  M-  and 

r D S 

tfg.  An  application  of  the  codes  DATCON  and  MAGIC  to  a large 
population  of  Eurasian  events  is  described  in  the  next  section 
of  this  report. 
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IV.  SPECTRAL  MAGNITUDE  DISCRIMINANT 
4.1  TECHNIQUE 

[71 

A variable  frequency  magnitude  discriminant1 
employing  spectral  magnitudes  m^  computed  from  short-period 
compress ional  waves  is  presently  being  tested  using  a large 
population  of  Eurasian  earthquakes  and  underground  explosions 
recorded  at  LASA^22^  and  the  original  Oyer  subarray  at  Norway. 
The  Norway  results  are  still  preliminary  and  will  be  dis- 
cussed in  a later  report. 

This  discriminant  is  designed  to  exploit  spectral 
differences  between  earthquakes  and  explosions.  Essentially, 
it  consists  of  a comparison  of  a spectral  magnitude  measure- 
ment in.  (f  ) at  a relatively  low  frequency  (e.g.,  f * 0.4  Hz) 

D l __  1 

to  a measurement  mb(f  ) at  a relatively  high  frequency  (e.g., 

f =2.5  Hz).  The  spectral  body  wave  magnitude  inb  is  defined 
as  the  log  of  the  amplitude  of  the  output  of  a very  narrow 
band,  phaseless  filter  (high  Q filter),  centered  at  frequency 
f plus  a distance  correction  factor  F.  That  is, 

5b  * 108  CW  * 

where  the  center  period  T^  = 1 / f g and  Ag  denotes  the  ampli- 
tude of  the  filter  output.  A particularly  desirable  feature 
of  this  type  of  discriminant  is  that  its  application  depends 
solely  on  the  recording  of  compressive  waves  from  earth- 
quakes and  explosions.  By  obviating  the  requirement  for  the 
recording  of  surface  waves , the  spectral  magnitude  technique 
should  be  applicable  to  smaller  magnitude  events  than  the 
conventional  mb  - Mg  discriminant. 
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4.2  EVENTS  RECORDED  AT  LASA 

The  event  population  employed  in  the  first  test  of  the 
spectral  magnitude  discriminant  consists  of  34  presumed  under- 
ground explosions  and  156  earthquakes  recorded  at  LASA.  Most 
of  these  events  occurred  in,  or  very  close  to,  the  USSR  or 
China,  at  epicentral  distances  of  60°  to  100°  from  LASA.  It 
is  important  to  note  that  only  about  55  percent  of  the  earth- 
quakes in  this  data  base  were  reported  in  the-  Monthly  Listings 
of  events  published  by  the  USCGS.  Thus,  45  percent  of  the 
earthquakes  in  this  data  base  are  of  either  unknown  or  poorly 
determined  focal  depth. 

Figures  4.1  and  4.2  are  plots  of  spectral  magnitudes 
at  several  different  frequencies  (low  frequency  m^  versus  high 
frequency  mb)  for  all  34  explosions  (closed  circles)  and  the 
shallow  earthquakes  (open  circles)  in  the  data  base.  Figure 
4.1  clearly  demonstrates  the  general  enrichment  of  the  high 
frequency  portion  of  the  explosion  body  wave  spectra  relative 
to  the  earthquake  spectra  over  the  entire  magnitude  range  of 
the  data.  For  instance,  for  a given  value  of  mb  (0.5  Hz)  the 
explosions  exhibit  mb  (2.0  Hz)  values  that  are  typically  0.6 
to  0.9  units  larger  than  the  earthquake  mb  (2.0  Hz)  values. 

In  particular,  note  that  this  separation,  or  discrimination, 
of  earthquakes  and  explosions  holds  over  the  entire  magnitude 
range  with  no  indication  of  convergence  of  the  two  populations 
at  small  magnitudes. 

Comparison  of  Figs.  4.1  and  4.2  shows  the  dependence 
of  the  spectral  magnitude  discriminant  on  parameters  such  as 
background  noise,  instrument  response,  and  the  frequency 
separation  of  the  magnitude  determinations.  First  note  the 
intermingling  of  event  populations  in  Fig.  4.2  at  small  magni- 
tudes, while  as  indicated  in  Fig.  4.1  only  one  small  magnitude 
(LASA  mb  = 3.9)  explosion  lies  in  the  earthquake  population. 
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Seismograms  for  this  explosion  (Fig.  4.3)  and  for  two  other 
explosions  with  in^  (0.5  Hz)  values  <4.5  clearly  demonstrate 
that  the  very  low  s ignal- to-noise  ratio  (S/N)  is  the  reason 
for  the  failure  to  discriminate  this  event.  Comparison  of 
the  seismogram  of  the  small  m^  = 3.9  explosion  (Fig.  4.3)  with 
the  seismogram  for  a larger  magnitude  (m^  ■ 5.8)  and,  there- 
fore, better  recorded  explosion  in  terms  of  S/N  (Fig.  4.4), 
indicates  that  the.  spectral  magnitude  determinations  based 
on  the  small  explosion  time  series  are  essentially  measure- 
ments of  the  noise  background,  in  particular  noise  in  the 
0.4  to  0.5  Hz  range. 

The  fact  that  many  of  the  small  magnitude  events  that 
fail  to  discriminate  in  Fig.  4.2  discriminate  well  in  Fig. 

4.1  is  a result  of  narrowing  the  frequency  spread  of  the 
magnitude  estimates  to  better  coincide  with  the  amplitude 
response  characteristics  of  the  LASA  instrumentation  and  the 
spectral  behavior  of  the  prevailing  background  noise.  This 
interpretation  is  supported  by  the  enhanced  discrimination  of 
larger  magnitude  events  in  Fig.  4.2  compared  to  Fig.  4.1,  in 
particular  for  earthquakes  6,  12  and  17.  In  this  case  S/N  is 
high  over  a wider  frequency  range  and  the  spectral  magnitudes 
computed  at  the  more  widely  separated  frequencies  yield  the 
enhanced  discrimination  observed  in  Fig.  4.2. 

It  has  been  pointed  out ^ that  maximum  discrimination 
of  event  populations  are  expected  for  spectral  magnitudes 
in^ ( f ) computed  at  frequencies  of  maximum  separation.  The 
actual  frequencies  will  depend  upon  microseismic  noise  spec- 
tra and  detector  site  characteristics.  In  addition,  in  order 
to  avoid  the  effects  of  spectral  maxima  and  minima  that  occur 
for  both  earthquakes  and  explosions,  it  is  probably  necessary 
to  use  combinations  of  spectral  magnitudes  at  several  different 
frequencies . 
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4.3  MULTIPLE  EVENTS 

The  variable  frequency  magnitude  discriminant  is  es- 
pecially suited  for  discrimination  of  multiple  explosion 
sequences  (Figs.  4.5  and  4.6)  that  are  designed  to  appear 
earthquake- like  in  terms  of  conventional  (m^  - Ms,  depth  of 

focus,  first  motion)  discriminants.  A multiple  event 

T 231 

scenario,  similar  to  one  previously  proposed,1  1 was  de- 
vised by  superposing  eight  seismograms  of  a presumed  explo- 
sion in  Eastern  Kazakhstan  recorded  at  LASA.  The  actual 
explosion  signature  used  in  this  scenario  is  shown  in  Fig. 
4.4.  The  array  of  explosions  and  their  relative  yields  were 
designed  to  produce  an  earthquake- like  seismogram  over  at 
least  some  range  of  azimuth.  The  particular  array  configura 
tion  (spacing  and  firing  order)  is  indicated  in  the  lower 
right  land  corner  of  Fig.  4.5.  Each  seismogram  comprising 
the  multiple  event  was  delayed  in  time  relative  to  the  first 
and  scaled  in  amplitude  to  correspond  to  an  explosion  of  the 
indicated  yield. 

Spectral  magnitude  estimates  of  the  multiple  event 
(closed  triangles)  are  plotted  in  Figs.  4.5  and  4.6  together 
with  the  event  populations  previously  plotted  in  Figs.  4.1 
and  4.2,  respectively.  The  numbers  1-5  indicate  the  dif- 
ferent azimuths  to  the  recording  sites  shown  in  Fig.  4.5. 
Note  that  the  legend  is  not  drawn  to  scale  as  the  shot 
spacing  is  1 km  and  the  recording  sites  are  at  teleseismic 
distances  ('v  83°) . 

The  multiple  event  sequence  was  especially  designed 
to  appear  earthquake-like  at  stations  broadside  to  the  shot 
array.  In  particular,  the  composite  seismogram  (Fig.  4.7) 
at  these  azimuths  exhibits  an  apparent  rarefaction  for  the 
first  prominent  motion.  It  should  be  noted  that  the  P-wave 
train  in  Fig.  4.7  is  not  preceded  by  background  noise.  If 
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Figure  4.5.  Spectral  magnitude  estimates  for  the  multiple 
event  sequence  and  the  Eurasion  event  population  plotted  in 
Fig.  4.1.  The  individual  shots  in  the  array  are  evenly 
spaced  and,  in  the  computations,  were  taken  to  be  1 km  apart 
The  shot  yields  are:  1 at  4 kT;  2 at  12  kT ; 3,  4 and  5 at 

31  kT;  6,  7 and  8 at  100  kT,  corresponding  to  amplitude  seal 
ing  of  1,  3,  7.8,  7.8,  7,8,  25,  25  and  25. 
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it  were,  the  very  small  excursion  starting  at  about  1 second 
would  most  likely  be  obscured  resulting  in  the  larger  rare- 
faction at  about  1.7  second  being  identified  as  the  first 
motion.  The  composite  seismogram  would  undoubtedly  lead  an 
analyst  to  measure  amplitudes  for  m^  based  on  one  of  the 
earlier  explosions  of  smaller  yield  while  measuring  amplitudes 
for  Mg  that  would  be  the  superposition  of  the  three  high 
yield  explosions.  The  net  effect  of  these  measurements  would 
be  to  increase  Ms  relative  to  mb  making  the  multiple  event 
approach  or  lie  within  the  earthquake  population  on  an 
mb  " Ms  Plot'  sPectral  magnitude  technique,  on  the 

other  hand,  by  basing  an  inb  measurement  on  the  largest 
amplitude  arrit-al  at  some  frequency  f in  the  wavetrain, 
exhibits  pronounced  discrimination  (Figs.  4.5  and  4.6)  of 
the  multiple  event  at  all  azimuths. 

The  composite  seismogram  in  Fig.  4.7  is  dominated 
by  the  three  largest  explosions  occurring  toward  the  end  of 
the  proposed  event  scenario.  In  order  to  simulate  a more 
earthquake-like  seismogram  a second  multiple  explosion 
scenario  was  devised.  Once  again,  the  presumed  explosion 
signature  shown  in  Fig.  4.4  was  used  as  the  primary  signal. 
Most  of  the  shot  array  characteristics  (spatial  separation, 
firing  order,  and  time  delays)  were  identical  to  the  pre- 
vious experiment,  the  major  difference  being  in  the  amplitude 
(yield)  ratios  of  successive  explosions.  In  this  second 
multiple  explosion  experiment  the  amplitude  ratios  were  1, 

3.1,  5,  10,  10,  20,  15  and  12.5.  The  primary  intent  of  this 
choice  scaling  factors  was  to  simulate  more  emergent 
composite  seismograms  devoid  of  suspicious  prominent  arriv- 
als. The  time  delays  and  scaling  were  still  such  that  Mg 
for  the  multiple  event  would  be  enhanced  relative  to  mb . 
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The  resulting  composite  seismograms  are  shown  in 
Fig.  4.8  at  five  different  azimuths  with  respect  to  the 
shot  array.  The  first  point  to  be  noted  is  that  with  the 
addition  of  noise  to  the  beginning  of  each  seismogram  the 
first  motions  at  azimuths  1,  3 and  3 would  most  probably 
be  picked  as  rarefactions.  Secondly,  the  complexity  of 
each  composite  signal  has  been  greatly  increased  over  that 
of  the  primary  signal  shown  in  Fig.  4.4.  An  analyst  making 
amplitude  measurements,  to  be  used  for  m^  determinations  in 
the  conventional  manner  (maximum  amplitude  within  the  first 
3 or  4 cycles),  would  undoubtedly  pick  amplitudes  correspond- 
ing to  the  earlier  smaller  explosions  in  the  sequence.  This 
is  equivalent  to  a further  enhancement  of  Mg  and  results  in 
an  event  that  appears  to  be  an  earthquake  on  an  Mg  - mb 
basis.  However,  when  subjected  to  the  spectral  magnitude 
technique  the  five  composite  signals  shown  in  Fig.  4.8  were 
found  to  discriminate  in  exactly  the  same  manner  as  the 
previous  scenario  (Figs.  4.5  and  4.6). 
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Figure  4.8.  Composite  seismograms  written  at  five  different 
azimuths  for  the  revised  multiple  event  scenario.  Approxi- 
mately 30  seconds  of  signal  is  shown  in  each  case;  the  two 
enlarged  tick  marks  along  each  abscissus  represent  a time 
increment  of  4.6  seconds.  The  maximum  O-P  amplitudes  of  each 
signal  are  about  50  mia. 


a mm 


KlfWOT  "W  V 

R-2483 


V.  CONCLUDING  REMARKS 

Most  of  the  important  conclusions  concerning  the 
variation  of  teleseismic  magnitude  with  the  various  contri- 
buting parameters  are  set  forth  in  Appendix  A and  will  not 
be  repeated  here.  However,  recent  work  has  revealed  the 
important  role  that  crustal  reverberations  may  play  in 
shaping  the  P wavetrain.  For  many,  though  not  all,  source 
regions  and  receiver  stations,  these  effects  should  be 
included  in  the  calculations. 

The  collection  and  playing  out  cf  seismic  data  from 
a variety  of  sources  has  been  a not  inconsiderable  task 
accomplished  during  the  first  portion  of  the  contract  period. 
Due  to  magnetic  tape  I/O  incompatibilities  between  various 
computer  systems,  the  initial  interaction  with  a data  source 
is  often  difficult.  However,  this  interaction  is  now 
routine  with  most  important  data  repositories. 

The  spectral  magnitude  discriminant  suggested  by 
Archambeau  has  been  extensively  tested  using  Eurasian  earth- 
quakes and  presumed  explosions  and  appears  to  offer  a 
distinct  improvement  over  previously  suggested  discrimination 
techniques  (e.g.,  m^  - Mg) . 
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APPENDIX  A 


EXPLOSION  YIELD  ESTIMATES  BASED  ON  SEISMIC  MAGNITUDE 
OBSERVATIONS:  DEPENDENCE  OF  MAGNITUDES  ON  YIELD, 

NEAR  SOURCE  MATERIAL  PROPERTIES,  DEPTH  OF  BURIAL 
AND  EARTH  STRUCTURE 


I.  INTRODUCTION 


We  consider,  from  the  onset,  that  the  simplest  and 
most  reasonably  implemented  approach  to  the  estimation  of 
explosive  yield  using  seismic  information  is  through  the 
direct  use  of  the  observed  seismic  magnitude  data  for  the 
event.  We  assume  that  seismic  magnitudes  will  be  measured 
from  a number  of  detectors  in  the  teleseismic  range  beyond 
1000  km  and  over  a range  of  azimuths  covering, at  least,  an 
interval  of  180°. 

The  magnitude  data  that  we  will  consider  are  m^,  or 
body  wave  magnitude,  and  Mg,  or  surface  wave  magnitude, 
observed  for  a given  event  at  a number  of  detectors.  In 
the  present  discussion  we  will  concentrate  most  of  our  at- 
tention on  the  body  wave  magnitude  data  since  the  P wave 
signal  from  which  mb  is  calculated  is  the  most  easily  detected 
at  teleseismic  distances  and  is  most  simply  interpreted  in 
terms  of  yield.  In  particular,  the  first  cycle  or  two  of 
the  wave  signal  is  less  perturbed  by  source  zone  "anisotropies" 
or  asymmetries  such  as  oriented  joints  in  the  rock,  initial 
prestress,  or  the  initial  cavity  shape.  Further,  earth 
structure  variations  between  the  different  source-receiver 
paths  are  strongest  in  the  crust  and  upper  mantle  and  there- 
fore will  usually  effect  surface  waves  more  strongly  then 
teleseismic  body  waves.  To  some  degree,  however,  the  stabi- 
lity of  the  mb  measurement  is  similarly  reduced,  since  the 
high  frequency  body  wave  is  subject  to  "scattering"  by  the 
crustal  layering  near  the  source  and  at  the  different  re- 
ceivers, giving  rise  to  fluctuations  from  receiver  to 
receiver.  Further,  the  anelastic  effects  are  much  stronger 
for  the  high  frequency  body  waves  than  they  are  for  the  much 
lower  frequency  surface  waves.  This  rather  strong  dependence 
of  observed  m^  values  on  structural  detail  at  the  source  and 
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receiver  and  on  the  upper  mantle  anelastic  properties  requires 
that  we  have  reasonably  detailed  knowledge  of  the  relevant 
eartli  structure.  However,  it  is  felt  that  corrections  for 
the  variations  of  structure  that  effect  the  body  wave  magni- 
tude measurements  are  more  easily  and  reliably  affected  tnan 
are  the  corrections  required  for  surface  waves.  Hence,  we 
have  concluded  that  an  mb  measurement,  followed  by  earth 
structure  corrections,  will  probably  give  the  best  yield  esti- 
mates. Even  so,  the  complementary  use  of  data  is  highly 
desirable,  since  it  provides  an  independent  measurement  of 
yield  that  can,  under  favorable  circumstances,  be  as  reliable 
or  even  more  reliable  than  the  body  wave  magnitude  data. 

In  particular,  in  a very  careful  study  of  seismic  magni- 
tudes of  underground  explosions,  Basham  and  Horner  (1973)  con- 
cluded that  Rayleigh  wave  magnitude  (Ms)  is  a more  consistant 
indicator  of  explosion  yield  than  P wave  magnitude.  The  M 
formula  used  by  Basham  and  Horner  is 

M * log  Ad  + B'  (A)  ♦ P(T)  ♦ S,  (l.D 

which  was  suggested  by  Marshall  and  Basham  (1972)  . In  this 
formula,  AR  is  the  maximum  observed  Rayleigh  wave  displace- 
ment, B ' (A)  is  the  distance  normalization  term  and  P(T)  is  a 
path  correction  term  which  is  a function  of  the  period  (T) 
at  which  AR  is  measured.  The  quantity  S is  a station  cor- 
rection determined  by  a statistical  analysis  of  many  explo- 
sion events  recorded  by  the  Canadian  Seismic  Network.  These 
values  appear  to  be  the  best  currently  available  in  terms 
of  consistancy  of  measurement  and  minimization  of  path  and 
station  effects. 

Basham  and  Horner  plotted  Mg  versus  yield  for  11  NTS 
explosions  in  tuff  and  rhyolite  along  with  two  (Gasbuggy  and 
Rulison)  explosions  in  western  Colorado.  The  Mg  values  for 
an  additional  set  of  18  NTS  tuff/rhyolite  explosions  with 


56 


unpublished  yields  were  also  gi'fen.  A least  squares  linear 
fit  for  the  explosions  with  published  yield  was  given  by 
Basham  and  Horner  as 

Mg  * 1.23  log  W + 1.56  . (1*2) 

Adding  the  yield  data  for  the  remaining  events,  this  relation 
is  altered  to 

* 1.19  log  W + 1.59  . (1*3) 

Calculating  the  error  by  dividing  the  difference  be- 
tween predicted  and  actual  yields  by  the  actual  yield,  one 
finds  that  Eq.  (1.3)  gives  an  average  error  of  30  percent 
for  these  28  events.  Several  of  the  events  have  anomalously 
large  errors.  Removing  the  four  worst  of  these  from  con- 
sideration, the  average  error  is  reduced  to  23  percent. 

Using  the  same  approach  for  body  wave  magnitude  obser- 
vations, Basham  and  Horner  established  the  empirical  magnitude 
yield  relationship 


mg  = 1.09  log  7 + 3.07  U 

where  mg  is  defined  as 

mg  = mb  + log  T ( 1 

with  T the  period  of  the  cycle  of  P wave  motion  used  to 
compute  the  magnitude  and  m^  the  usual  body  wave  magnitude 

defined  by 

mb  55  log  (A/T)  + B (1 
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Here  A is  the  zero  to  peak  ground  motion  in  millimicrons 
measured  from  the  largest  cycle  of  motion  in  the  first  few 
cycles  (i.e.,  first  three  cycles)  of  P wave  motion  fro* 
the  event.  The  B factor  is  the  distance  correction  based  on 
the  source  to  receiver  distance. 

Adding  yield  data  unavailable  to  Basham  and  Horner, 
the  relation  (1.4)  is  altered  to 

m*  = 0.913  log  7 + 3.53  . (1*’ 

D 

Once  again  the  error  may  be  calculated  by  dividing  the  dif- 
ference between  predicted  and  actual  yields  by  the  actual 
yield.  Thus,  one  finds  that  Eq.  (1.7)  gives  an  average 
error  of  50  percent  for  these  28  events*.  This  should  be  com 
pared  to  the  30  percent  average  error  in  yield  prediction 

based  on  Ms« 

f is  of  considerable  importance  that  we  isolate  the 

origins  of  the  scatter  in  the  magnitude-yield  relation  since 
even  for  the  carefully  measured  and  empirically  corrected 
data  from  Basham  and  Horner,  a rather  large  scatter  remains. 

I This  results,  of  course,  in  a rather  large  uncertainty  in 

estimation  of  yield  if  these  purely  empirical  relations 
are  used.  It  is  clear,  however,  that  we  will  have  different 
magnitude-yield  relations  for  explosions  in  different  near 
r source  materials  and  if  the  receiver  array  is  such  as  to  sample 

the  seismic  field  after  it  has  propagated  through  vastly 
different  earth  structures.  Thus,  the  magnitude-yield 
relation  for  hard  rock  materials  is  different  than  that  for 
t soft  rock  materials  and  further  the  relation  will  be 

different  as  measured  by  an  array  in  a shield  area  as 
opposed  to  a detector  array  in  an  active  tectonic  area. 

It  is  desirable,  therefore,  to  establish,  in  a deter- 
* ministic  fashion,  magnitude-yield  relations  that  are  appro- 

priate for  particular  receiver  arrays  and  source  locations 
where  depth  and  near  source  material  properties  are  also 
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variable.  Equivalently,  the  magnitude  observations  may  be 
corrected  for  source-receiver  structure  effects,  near  source 
material  properties,  depth  and  other  variables  affecting  the 
magnitude  observations  in  a deterministic  fashion  in  order 
to  greatly  reduce  the  scatter  in  a magnitude-yield  relation. 

In  either  approach  it  is  necessary  to  determine  which  para- 
meters affect  the  magnitude  measurement,  and  the  precise  dependence 
of  the  magnitude  measurement  on  these  parameters.  Therefore 
in  this  report  we  will  focus  our  discussion  on  the  variability 
of  mb  and  Mg  with  respect  to  such  variables  as  near  source 
material  properties,  initial  medium  stress  state,  earth 
structure  and  source  depth. 

In  order  to  see,  in  a simple  and  direct  way,  the 
relationship  of  seismic  magnitude  to  the  source  properties, 
we  will  introduce  the  reduced  displacement  potential,  ^(t). 

We  can  easily  demonstrate  that  all  seismic  magnitudes  depend 
on  the  single  source  function,  \p,  in  a simple  way. 

Consider  first  the  seismic  body  wave  magnitude  mb  from 
an  explosion.  For  a spherically  symmetric  explosion  the  reduced 
displacement  potential,  'J'(t),  is  related  to  the  radial  displace- 
ment, u (t,r),  by 


u(t.r)  - i&l  . 


(1.8) 


where  t is  time, 


r is  radial  distance  from  the  source, 


a is  compressional  velocity  in  the  source  region, 


x ■ t-r/a,  the  retarded  time. 


Note  that  the  steady  state  value  of  the  reduced  displacement 
potential  (RDP)  is  given  by 


’M00)  3 r2u(“,r) 


X 


where  u(“,r)  is  the  static  displacement  at  the  radius  r. 

Since  the  seismometer  filters  the  displacement  wave 
form,  it  is  convenient  to  transform  Eq.  (1.8)  into  the 
frequency  domain.  We  find 


u 


(1.10) 


In  the  teleseismic  field,  many  wavelengths  from  the  source, 

2 

the  term  attenuating  as  r"  is  negligible.  Therefore,  the 
teleseismic  body  wave,  u^,  is 
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(1.11) 


Since  the  spectrum  of  ip  for  an  explosion  (see  Section  II) 
is  essentially  flat  in  the  teleseismic  frequency  band 
(f  < 2 Hz),  then  Eq.  (1.11)  may  be  written 


(1.12) 


where  we  have  used  a b sic  property  of  the  Fourier  transform, 
namely 

A 

<H0)  * *(»)  • (1.13) 


Since  i K00)  bas  dimensions  of  L3  (Eq.  (1.9)),  then 
using  cube-root  scaling  we  find  that  both  iK00)  and  tele- 
seismic body  wave  displacement  should  scale  directly  as  the 
device  yield.  This  important  result  assumes  that  the  near 
source  material  properties  do  not  change,  the  only  variable 
being  the  yield  of  the  explosive. 


with 


Consequently  we  have  that 


mb  * log 
*(co)^7,  7 


being  the  explosion  yield 


(1.14) 
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For  the  surface  wave  magnitude  Mg  we  have  a simil 
result.  In  particular,  the  spectrum  of  the  vertical  s ~ace 
displacement,  # for  the  Rayleigh  waves  propagating  in  a 
multilayered  medium,  excited  by  a source  of  pure  dilatation 
in  a layer  with  material  index  s is  given  by  Harkrider 
(1964,  1971,  as 

wq  * - 4irius  4>(w)  KR  AR  eiTs  (kRr)  (1.15) 

where 

y = rigidity  of  the  material  in  layer  s, 

ip  (w)  = Fourier  transformed  reduced  displacement 
potential , 

kR  = u/cR  where  cR  is  the  Rayleigh  wave  velocity 
at  the  frequency  to, 

Kd  = depth  dependent  Rayleigh  wave  excitation 
R 

factor, 

Ad  * amplitude  response  of  the  multilayered  medium, 

T - UR  /a  where  R is  the  effective  "cavity" 

S 0 S 0 . 

radius  for  the  explosion  and  is  the  P 

wave  velocity  in  layer  s, 

1 1 C 2 ) (k  r)  = Spherical  Hankel  function  of  the  second  kind, 
o R 

Except  for  the  small  shift  in  origin  time,  tg,  the 
only  source  dependent  term  is  i|>s(w)»  the  reduced  displace- 
ment potential.  The  remaining  terms  are  nearly  independent 
of  depth  of  burial  for  depths  typical  of  underground  explo- 
sions. Further,  only  a small  correction  is  required  for  the 
case  of  an  explosion  in  one  medium  (e.g.,  tuff)  exciting 
Rayleigh  waves  in  an  average  source-receiver  crustal  struc- 
ture which  is  predominantly  of  gran:te-like  materials. 
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For  the  long  periods  which  are  of  interest  for  sur- 
face wave  observations,  the  source  spectrum,  4<s(w),  is 
equal  to  the  constant  value  (|>s(“)  (see  Section  II  also). 
With  these  observations,  (1.15)  implies  that  surface  wave 
magnitude  is  simply  related  to  the  source  by  the  pro- 
portionality relation 

M * log  [u  iH")  ] • ^ 

s 


We  see  that  the  magnitude  mb  i nd  Mg  should  therefore  be  re 
lated  to  yield  by  relations  of  the  form 


n»b  * log  T ♦ B 1 


(1.17) 


Ms  ^ l°g  I ♦ C ) 

where  the  coefficient  of  the  log  7 term  should  be  close 
to  unity  and  B and  C will  be  explicitly  dependent  on  the 
compressional  velocity  and  rigidity  of  the  medium  at  the 
source,  respectively. 


Therefore,  we  see  that  all  source  parameter  effects 
on  seismic  magnitudes  may  be  discussed  in  terms  of  iKu), 

\x  and  a . 


The  success  of  a deterministic  approach  to  magnitude 
prediction,  wherein  earth  structure,  near  source  material 
properties  and  a variety  of  other  factors  are  taken  into 
account  can  be  gauged  from  the  results  shown  in  Figure  1.1. 
Here  the  circles  indicate  m^1-*  - Mg  values  computed  from 
theoretically  produced  seismograms,  where  mb  is  the 
magnitude  measured  from  the  first  cycle  of  F wave  motion  and 
M is  measured  from  the  20  second  Rayleigh  wave  (vertical  com- 
ponent). The  line  to  the  right  of  the  circles  indicates  the 
mb  value  obtained  from  a measurement  of  the  largest  P wave 
amplitude  in  the  wave  train  (the  second  cycle  of  motion  was 
actually  used).  The  squares  indicate  the  values  of  Mg  that 
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Figure  1.1.  Rayleigh  surface  wave  magnitude  as  a function 
of  body  wave  magnitude  for  explosions  in  granite  (model  G-l). 
The  circles  denote  theoretical  values  of  m^l)  and  M§  for  zero 
prestress,  the  squares  for  165  bar  prestress.  The  horizontal 
lines  indicate  the  range  of  possible  mb  values,  the  vertical 
lines  the  range  of  possible  M§  values.  The  upper  and  lower 
bound  lines  for  earthquake  data  are  shown  for  comparison  along 
with  the  mean  NTS  explosion  magnitude  line.  The  X symbols 
are  individual  explosions  defining  the  extreme  upper  range 
of  observed  Ni§  values  for  explosions.  The  yields  of  the  ex- 
plosions are  indicated  along  the  m^,  scale. 
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are  measured  for  the  same  explosions  but  with  tectonic 
release,  under  the  condition  that  the  medium  was  prestressed 
in  shear  to  a 165-bar  level.  The  vertical  lines  indicate 
the  range  of  Mg  that  could  be  assigned  to  the  individual 
explosions,  the  precise  value  depending  on  the  period  of  the 
surface  wave  used  in  the  measurement.  The  extreme  (upper) 
value  of  the  Mg  magnitude  for  a given  explosion  is  obtained 
if  the,  so  called,  Airy  phase  of  the  surface  wave  train  is 
used,  corresponding  to  an  amplitude  measurement  at  a period  of 
from  12  to  14  seconds.  Clearly  then,  measurement  procedures 
will  affect  the  magnitudes  obtained  and  a stricter  definition 
of  the  magnitude  measurement  is  called  for  in  order  to  reduce 
the  scatter  in  m^  and  Mg  measurements.  (No  latitude  in  the 
cycle  used  for  the  P wave  amplitude  measurement  and  only  the 
20-second  surface  wave  measured,  for  example.) 

These  results  are  in  good  agreement  with  observations, 
as  is  shown  by  the  NTS  » mb  - Mg  line.  The  crosses  appearing 
in  the  figure  are  large  Mg  explosions  and  indicate  the  upper 
bound  on  the  m'^-Mg  population.  The  fact  that  the  theoretical 
results  are  somewhat  above  the  NTS  line  is  expected,  in  that 
the  medium  model  used  was  for  granite  while  most  of  the  NTS 
explosions  were  in  soft  materials  (e.g.,  tuff). 

The  results  may  be  put  in  the  form  of  theoretically 
predicted  magnitude-yield  relations  which  are  shown  in  Figures 
1.2  and  1.3.  These  may  be  compared  with  the  empirically 
determined  equations  (1.3)  and  (1.7).  The  agreement  is 
quite  good  considering  the  fact  that  the  empirical  data 
includes  explosions  in  a variety  of  near  source  material 
types . 
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These  theoretical  results  strongly  suggest  that  the 
detailed  source  description  used  is  accurate  and  that  we  can 
use  this  capability  to  either  predict  the  magnitude-yield 
relations  for  a wide  variety  of  source  and  medium  conditions 
or  that  we  can  correct  observed  magnitude  data  for  variability 
dther  than  yield,  to  obtain  very  accurate  yield  estimates. 

In  the  following  sections  we  give  a detailed  discussion 
of  all  the  factors  influencing  magnitude,  from  the  point  of 
view  of  yield  estimation  based  on  magnitude,  with  each  of  these 
sections  containing  a concise  summary. 
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II.  THE  DEPENDENCE  Or  TELESEISMIC  MAGNITUDES  ON 
NEAR  SOURCE  MATERIAL  PROPERTIES 

2.1  INTRODUCTION 

In  order  to  determine  teleseismic  amplitude  dependence 
on  near  source  rock  environment,  it  is  necessary  to  obtain  an 
equivalent  elastic  source  that  matches  the  propagation  charac- 
teristics of  the  real  source  at  distances  over  which  the 
medium  response  is  linear  elastic.  Explosion  source  function 
prediction  is  possible  if  one  is  able  to  simulate  a propagat- 
ing stress  wave  in  at  least  a spherically  symmetric  geometry. 
For  nuclear  sources,  this  simulation  requires  an  understanding 
of  material  behavior  from  a few  megabars  in  pressure  down  to 
the  elastic  level.  The  source  function  is  obtained  from  the 
calculated  displacement  or  velocity  history  of  a particle  in 
the  elastic  region. 

By  and  large,  seismologists  have  been  forced  to  use 
experimentally  determined  explosion  source  functions  for 
teleseismic  ground  motion  predictions.  These  source  func- 
tions are  limited  in  number,  quality  and  the  rock  environment 
they  represent.  These  limitations  are  due  mainly  to  the  ex- 
pense and  technical  difficulty  involved  in  fielding  an  experi- 
ment designed  to  obtain  free  field,  time  history  data  in  the 
elastic  regime  for  a nuclear  source. 

Computational  techniques  and  material  models  have  been 
developed  at  S3  which  permit  a prediction  of  the  equivalent 
elastic  source  given  the  yield  of  the  explosive  and  the 
material  properties  of  the  near  source  geologic  environment. 

Since  laboratory  tests  on  rock  samples  provide  the 
essential  data  base  for  a preshot  estimate  of  material  proper- 
ties, it  is  important  to  determine  the  sensitivity  of  the 
equivalent  source  to  those  material  properties  which  are 


68 


i 


■ 


subject  to  routine  laboratory  measurement.  In  this  chapter 
we  present  the  results  of  a computational  parameter  study 
which  identifies  the  dependence  of  teleseismic  magnitude  on 
near  source  material  properties. 

2.2  MODEL  CONFIRMATION  AND  RELEVANCE  TO  TELESEISMIC  GROUND 
MOTION  

Befcre  presenting  the  results  of  the  parameter  study, 
it  xs  important  first  to  establish  at  least  some  confidence 
in  the  material  models  capability  to  predict  close  in,  free 
field  ground  motion  and  second, to  identify  the  manner  in 
which  changes  in  material  parameters  will  be  related  to 
changes  in  teleseismic  ground  motion. 

Figure  2.1  shows  measured  and  predicted  ground  motion 
from  a nuclear  explosion  (CYPRESS)  at  ranges  of  550,  700, 
and  1100  feet  from  the  detonation  point.  Figure  2.2  compares 
measured  and  predicted  ground  motion  from  a chemical  explo- 
sion (Mine  Dust  HE)  at  12  feet  from  the  center  of  the  explo- 
sive. Figure  2.3  shows  measured  and  predicted  stress  history 
from  a nuclear  explosion  (Misty  North)  at  180  feet  from  the 
detonation  point.  While  agreement  is  certainly  less  than 
perfect  in  all  cases,  the  computational  technique  does  seem 
capable  of  predicting  peak  positive  and  negative  velocities 
along  with  the  overall  shape  of  the  propagating  stress  pulse. 
The  late  time  disagreement  between  calculated  and  measured 
ground  motion  for  CYPRESS  could  be  due  to  reflected  waves 
obscuring  the  free  field  portion  of  the  measurements. 

The  results  presented  in  Figures  2.1  through  2.3  show 
velocity  or  stress  histories  at  various  distances  from  the 
source.  For  teleseismic  ground  motion  it  is  possible  to  re- 
duce the  coupling  efficiency  to  a single  number,  ij^,  which 
scales  directly  with  the  yield  of  the  explosion. 
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Figure  2.3.  Stress-time  waveform  at  180  feet  and  predicted 
waveform  at  same  location  for  Misty  North. 
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For  a spherically  symmetric  source,  the  reduced  dis- 
placement potential,  ip,  is  defined  by 


u(t,r) 


(2.1) 


where  u(t,r)  * radial  particle  displacement  at  t and  r, 
a = P wave  velocity  in  the  medium, 
r * radial  position, 
t = t - r/a,  retarded  time, 

^ = dip/dT , reduced,  velocity  potential. 

In  the  elastic  region 


i !l  « 02  ill 

3t 2 3r2 

Therefore  * depends  only  on  t and  its  wave  form  does  not 
change  with  changing  r. 

Given  the  displacement  time  history,  u(t,r),  at  an 
aribtrary  radius  from  the  source,  then  Eq.  (2.1)  may  be  used 
to  determine  i|>(t)  . The  value  of  <Kt)  so  determined  will 
be  independent  of  radius  if  the  displacement  time  history 
was  chosen  in  the  elastic  zone. 

From  Eq.  (2.1)  the  steady  state  value  of  the  reduced 
displacement  potential  is  given  by 

(«=)  = r2u(“,r)  (2.2) 

where  u(°°,r)  is  the  static  displacement  at  the  radius,  r. 

For  all  the  calculated  reduced  displacement  potentials 
obtained  to  date  it  has  been  found  that  teleseismic  body  wave 
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magnitude,  mb>  has  the  proportionality 


m^  ^ log 


(2.3) 


while  surface  wave  magnitude,  Mg,  obeys 

M 'v  log  [y  C00) ] (2*4) 

s 

where  y is  the  shear  modulus  of  the  near  source  material. 

Equations  (2.3)  and  (2.4)  are  valid  for  device  yields 
up  to  approximately  2S0  kT.  Their  form  is  a consequence  of 
the  source  spectrum's  frequency  content  in  the  teleseismic 
band  (f  < 2 Hz).  Our  calculated  spectra  are  essentially  flat 
over  this  frequency  interval  for  device  yields  up  to  250  kT. 
Therefore, changes  in  teleseismic  magnitudes  caused  by  changes 
in  near  source  material  properties  may  be  obtained  by  com- 
paring the  values  of  tj>(»)  appropriate  for  each  material. 

If  i|).(«.j  and  <jJk(“)  are  two  values  of  tH00)  corres- 
ponding to  materials  i and  k then  the  change  in  tele- 
seismic magnitude  is  given  by 


Am  = m m^=  log 


^(“0" 
^k(°°)  ’ 


(2.5) 


where 


ak 

a = — for  m, 
a ■ b 


and 


a 


yi 

— for  M . 
^k 
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If  material  properties  remain  constant  and  the  only  change 
is  device  yield,  , then 


^kC")  = 57  C°°D  . 

and  Eq.  (2.5)  becomes 


Equation  (2.5)  gives  the  "scaling"  law  for  teleseismic 
magnitudes  in  terms  of  4(°°)-  This  equation  is  always  valid 
as  long  as  the  spectrum  of  the  equivalent  source  is  flat 
within  the  teleseismic  frequency  band  and  path  effects  asso- 
ciated with  the  earth  structure  are  invariant  between  events. 

2.3  THE  SENSITIVITY  OF  iH00)  TO  AIR-FILLED  POROSITY,  MATERIAL 
STRENGTH  AND  OVERBURDEN  PRESSURE 

A number  of  calculations  have  been  performed  in  order 
to  identify  the  sensitivity  of  4(°°)  to  changes  in  air-filled 
porosity,  material  strength  and  overburden  pressure.  Table 
2.1  lists  the  material  properties  along  with  the  calculated 
value  of  iK00)  for  a source  yield  of  0.02  kT . Also  included 
in  the  table  is  Am  computed  from  Eq.  (2.5)  where  calculation 
5 has  been  chosen  for  the  reference  magnitude,  i.e., 
k * 5 m Eq.  (2.5).  The  P wave  velocity  is  2.4  km/sec 
and  the  shear  modulus  is  40  kbars  for  all  calculations. 
Therefore,  the  value  of  "a"  in  Eq.  (2.5)  is  unity  for  this 
study. 

The  volume  fraction  of  air-filled  voids  in  the  material 
is  obtained  from  the  column  labeled  The  effect  of  air- 

filled  porosity  on  magnitude  is  obtained  by  comparing  the 
results  of  calculations  1,  6 and  4.  Figure  2.4  shows  Am 
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TABLE  2.1 

SUMMARY  OF  SOURCE  COMPUTATIONS 


n 


versus  for  these  three  calculstions . We  see  that  Am 

o 

is  highly  sensitive  to  changes  in  air-filled  porosity.  Over 

the  interval  0 < 4>  < 0 . 02 .magnitude  decreases  by  0.313  units. 

o 

A material  with  2 percent  porosity  will  give  the  same  magni- 
tude as  a material  with  no  porosity  if  the  device  yield  is 
increased  by  i factor  of  two  over  that  used  at  the  zero  poro- 
sity site. 

Figure  2.5  shows  the  influence  that  the  maximum  allowable 
material  strength,  Y , has  on  magnitude.  This  figure  compares 
the  results  of  calculations  4 and  5.  Increasing  Ym  by  a 
factor  of  two  causes  magnitude  to  decrease  by  0.192  magnitude 
units.  Therefore,  the  device  yield  could  be  increased  by  a 
factor  of  1.55  in  the  strong  material. 

These  calculations  show  that  a porous,  high  strength 
material  optimizes  teleseismic  decoupling.  Increasing  both 
the  porosity  and  material  strength  may  obscure  the  device 
yield  by  at  least  a factor  of  three  when  teleseismic  magni- 
tude is  used  to  infer  yield.  We  should  note  that  both  high 
air-filled  porosity  and  high  material  strength  imply  that 
the  material  characterizing  the  near  source  environment  is 
above  the  water  table. 

For  a given  value  of  Y , the  shape  of  the  material 
strength  surface  is  controlled  by  Pm,  the  pressure  at  which 
Y is  attained.  From  calculations  9 and  7 we  see  that  in- 

m 

creasing  Pm  by  a factor  of  2.08  causes  magnitude  to  in- 
crease by  0.26  units.  The  coupling  efficiency  increases 
as  P is  increased,  which  is  opposite  to  the  effect  of 
increasing  4>  or  Y . This  simply  indicates  that  the  shape 
of  the  failure  surface  plays  an  important  role  in  determining 
coupling  efficiency.  The  material  has  less  strength  at  low 
pressures  when  P^  is  increased  and  we  should  expect  a cor- 
responding incre.se  in  magnitude. 
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The  effect  of  overburden  pressure,  P , on  magnitude 
is  shown  in  Figure  2.6.  Coupling  efficiency  decreases  with 
increasing  overburden  pressure,  when  all  other  material 
properties  remain  constant.  An  increase  in  the  depth  of 
burial  (DOB)  by  a factor  of  six  causes  magnitude  to  decrease 
by  0.3  units,  which  in  turn  permits  a yield  increase  by  a 
factor  of  two.  The  decrease  in  magnitude  is  not  uniform 
with  depth.  An  increase  in  DOB  at  shallow  depths  causes 
magnitude  to  decrease  more  then  that  caused  by  an  equivalent 
increase  in  DOB  at  greater  depths.  This  effect  results  from 
the  rapid  increase  in  material  strength  at  low  pressures. 

We  should  note  that  it  is  unlikely  that  material 
properties  would  remain  invariant  as  overburden  pressure  is 
increased  from  100  kbars  to  600  kbars.  An  increased  DOB 
usually  corresponds  to  an  increase  in  water  content,  result- 
ing in  decreased  air-filled  porosity  and  material  strength. 
These  changes  would  more  than  compensate  for  the  overburden 
pressure  decoupling  effect. 

Calculations  9 and  10  show  that  magnitude  is  not  very 

sensitive  to  large  changes  in  Pc,  the  pressure  at  which  all 

ambient  porosity  is  irreversibly  removed.  An  increase  in 

p by  a factor  of  2.5  produces  an  increase  in  magnitude  of 
c 

0.03  units. 

Finally,  calculations  4 and  3 show  the  sensitivity  of 
magnitude  to  Pe,  the  pressure  limit  at  which  irreversible 
pore  collapse  begins.  An  increase  in  Pg  produces  an  in- 
crease in  coupling  efficiency.  A preshot  estimate  of  Pe 
that  is  incorrect  by  a factor  of  two  will  result  in  an  error 
in  predicted  magnitude  by  0,2  units. 
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2 . 4 SUMMARY 


Increases  in  air-filled  porosity,  material  strength, 
and  overburden  pressure  reduce  seismic  coupling.  Increas- 
ing the  pressure  at  which  irreversible  pore  collapse  begins 
produces  an  increase  in  seismic  coupling.  These  are  the 
critical  parameters  which  must  be  determined  prior  to  a con- 
fident prediction  of  teleseismic  magnitude.  If  any  of  these 
parameters  is  specified  without  appropriate  laboratory  tests 
being  performed  on  rock  samples  from  the  near  source  environ- 
ment, then  predicted  magnitude  should  be  suspect  b>  at  least 
0.2  magnitude  units. 

For  an  intrusive  igneous  rock,  i.e.,  granite,  we  would 
expect  the  major  material  properties  variation  to  occur  from 
changes  in  material  strength,  i.e.,  variations  in  porosity 
are  negligible.  Therefore,  the  variation  in  magnitude 
between  two  granite  sites  could  be  0.45  units  if  the  material 
strength  parameters  changed  by  a factor  of  two  between  the 
sites . 

The  variation  in  magnitude  from  changes  in  the  air- 
filled  porosity  parameters  is  0.513  units  where  we  have 
limited  the  porosity  to  be  less  than  2 percent  and  allowed 
for  a factor  of  two  variation  in  pressure  at  which  irreversible 
pore  collapse  begins.  For  sedimentary  and  volcanic  rocks 
both  porosity  and  material  strength  changes  between  sites 
could  produce  a corresponding  variation  in  magnitude  of 
approximately  1.0  magnitude  units. 

There  is  no  substitute  for  good,  preshot  material 
property  data  in  order  to  predict  the  effect  of  near  source 
geology  on  teleseismic  magnitude.  Laboratory  tests  on  ap- 
propriate rock  samples  coupled  with  density  and  P wave 
velocity  logs  at  the  site  are  able  to  furnish  the  necessary 
data  required  to  specify  the  material  properties  that 
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^ determine  seismic  coupling.  This  data  base  must  be  estab- 

lished if  the  effect  of  near  source  geology  on  teleseismic 
magnitude  is  to  be  predicted. 
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III.  DEPENDENCE  OF  SEISMIC  MAGNITUDES  ON  SOURCE 
DEPTH,  MEDIUM  PRESTRESS  AND  SOURCE  ASYMMETRY 

In  the  previous  section  the  influence  of  the  source 
material  parameters  on  the  reduced  displacement  potential,  iK°°)  » 
was  discussed.  This  quantity,  <K“)  * essentially  provides  the 
source  coupling  for  elastic  waves  emanating  from  the  source 
region  and  is  of  paramount  importance  for  determining  the 
amplitude  of  teleseismic  body  and  surface  waves.  However, 
other  features  also  influence  the  wave  train  from  which 
measurements  of  mb  and  Mg  are  ultimately  made.  One  is, 
of  course,  the  earth  struct  re  between  source  and  receiver 
which  is  the  subject  of  a subsequent  section  (Section  IV). 

In  this  section  we  consider  the  effects  of  medium  prestress, 
source  asymmetry  (e.g.,  cavity  shape)  and  free  surface 
generated  interference  effects,  which  are  dependent  on 

source  depth. 

3.1  THE  DEPENDENCE  OF  mb  AND  Mg  ON  DEPTH  OF  BURIAL 

The  burial  depth  of  the  explosion  source  affects  the 
teleseismic  P wave  record  in  several  ways.  By  far  the 
most  important  is  the  interference  between  the  direct  (P) 
and  free  surface  reflected  (pP)  waves.  Other  depth  of  burial 
(DOB)  effects  include  the  variation  of  local  material  proper- 
ties and  overburden  pressure.  A minor  consideration  for 
explosions  is  the  fact  that  large  DOB  changes  alter  the  source- 
to-receiver  travel  path  and  therefore  the  transfer  function. 

For  some  events  large  S waves  may  be  generated  due  to  tectonic 
release  or  other  causes.  In  this  case  the  DOB  related  time 
lag  between  P and  the  reflected  sP  phase  can  be  important. 
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of  course,  the  earth  structure  between  source  and  receiver 
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In  this  section  we  consider  the  effects  of  medium  prestress, 
source  asymmetry  (e.g.,  cavity  shape)  ana  free  surface 
generated  interference  effects,  which  are  dependent  on 

source  depth. 

3.1  THE  DEPENDENCE  OF  mb  AND  Mg  ON  DEPTH  OF  BURIAL 

The  burial  depth  of  the  explosion  source  affects  the 
teleseismic  P wave  record  in  several  ways.  By  far  the 
most  important  is  the  interference  between  the  direct  (P) 
and  free  surface  reflected  (pP)  waves.  Other  depth  of  burial 
(DOB)  effects  include  the  variation  of  local  material  proper- 
ties and  overburden  pressure.  A minor  consideration  foi 
explosions  is  the  fact  that  large  DOB  changes  alter  the  source- 
to-receiver  travel  path  and  therefore  the  transfer  function. 

For  some  events  large  S waves  may  be  generated  due  to  tectonic 
release  or  other  causes.  In  this  case  the  DOB  related  time 
lag  between  P and  the  reflected  sP  phase  can  be  important. 
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The  variation  with  depth  of  the  local  material  proper- 
tie  is  properly  considered  a source  effect  rather  than  a DOB 
effect.  The  dependence  of  the  source  spectrum  cn  local 
material  properties  was  discussed  in  Section  II. 

The  important  parameter  controlling  the  interference 
between  P and  pP  is  the  average  source- tc-surface  P wave 

cc  , -ifSS 

velocity  (V  ) . For  example,  a 10  percent  decrease  in  Vp 
results  in  an  11  percent  increase  in  the  arrival  time  lag  be- 
tween P and  pP  having  the  effect  of  increased  depth  of 
burial . 

As  the  overburden  pressure  (P^)  increases  with  depth, 
the  effect  is  to  decrease  the  average  low  frequency  spectral 
level  of  the  displacement  field  from  the  source.  For  example, 
in  Figure  3.1,  the  amplitudes  of  the  RVP  transforms,  |$»(w)|, 
from  three  SKIPPER  calculations  are  compared.  The  three  source 
models,  T-4,  T-5,  T-6  are  identical,  except  for  the  overburden 
pressure.  We  see  that  with  increasing  DOB,  the  RVP  transform 
exhibits  the  behavior  predicted  by  the  standard  scaling  laws; 
that  is,  a reduced  $(0)  and  a spectral  shift  to  higher 
frequencies.  However,  these  effects  are  quite  minor,  as  a 
more  than  two-fold  increase  in  depth  results  in  an  apparent 
yield  decrease  of  about  15  percent,  based  on  the  reduced  level 
of  4>(0).  (Note  that  ^(0)  = ip(°°).) 


3.2  INTERFERENCE  BETWEEN  P AND  pP 

The  effect  of  the  free  surface  on  the  spectrum  of  the 
observed  first  arrival  P wave  train  may  be  demonstrated  as 
follows . 

Let  u and  u „ denote  the  displacements  for  P 
P P*3 

and  pP,  respectively.  If 


Up  - F(t)  , 


(3.1) 


i 


_J 
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Figure  3.1.  Reduced  velocity  potential  transform 
amplitude  versus  frequency  for  tuff 
source  functions  T-4,  T-5  and  T-6. 

Scale  is  appropriate  to  a yield  of  0.02 
kt.  The  1 Hz  points  for  yields  of 
100  and  1000  kt  are  shown  for  reference. 
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then 


upP  = ‘AF(t‘T0)  * 


(3.2) 


where 


A = free  surface  reflection  coefficient  (%  1) , 

T % 2H/VSS  , 
o P 

H = depth  of  burial  (DOB) , 

C C 

V = average  source-surface  P wave  velocity. 


In  the  frequency  domain 


UP  + upP 


f(w)  ( 


-iwt 


1 - Ae 


0* 


(3.3) 


and 


Iup  * upp i 


• |F(«) I 


» |F(u) | 


1 - 2A  coswt  + A2 
o 


1 - A | , u>  H 2mr/T 


1/2 


t 

( I 1+A | , 


(3.4) 


(3.5) 


w = (2n+l)Tt/T  , n=  0,1,2,  ... 
0 


The  modulation  of  the  frequency  spectrum  is  given  by 
(3.5)  with  the  spectral  peaks  and  troughs  occurring  alternately 
at  frequencies 


nV 


SS 


f = — ^ , • • • • 


(3.6) 


87 


I—WNB  ■llli  ■■  '• 


The  effect  of  (3.6)  is  shown  in  Figure  3.2  where  |up|, 


u 


depths . 


and 


lu  + u r,  I are  plotted  for  two  yields  and  burial 
1 p pP  1 

The  two  source  functions  used  were  T-4  and  T-6  of 
Figure  3.1.  Thus,  with  the  exception  of  overburden  pressure, 


the  material  properties  in  the  source  region  were  held  constant 

ss 

for  both  calculations  and  V =2.4  km/sec. 


To  demonstrate  that  the  shift  of  the  spectral  peak 
occurring  for  the  examples  in  Figure  3.2  is  almost  entirely 
a DOB  effect,  we  show  in  Figure  3.3  the  vertical  component 
of  | Up  + Upp|  plotted  for  six  yields  ranging  from  800  kT 
with  the  source  to  surface  travel  time  held  constant. 


In  order  to  determine  the  effect  of  the  interference 

of  P and  pP  on  mb,  it  is  necessary  to  consider  the  time 

domain  signal  produced  by  the  two  phases.  In  Figure  3.4  the 

basic  wave  form  of  the  first  motion  P wave  train  for  shallow 

explosions  is  shown.  The  direct  P and  reflected  pP  waves  are 

shown  along  with  their  summation  as  recorded  by  a standard 

LRSM  short  period  seismometer.  In  Figure  3.5  the  sum  of 

P + pP  is  shown  for  a series  of  burial  depths.  For  all 

seismograms  of  Figures  3.4  and  3.5  the  source  T-3  and  the 

earth  structure  model  CIT  109  0.5Q  were  used.  The  epicentral 

distance  was  4087  km  corresponding  to  a ray  with  takeoff  angle 

SS 

0 = 13.8°.  For  each  of  the  calculations  Vp  = 2.4  km/sec 
and  the  DOB  was  varied.  However,  the  significant  parameter 
is  the  travel  time  lag  between  P and  pP, 


t.  " -ct C1  + cos  2e)  • 

L cose 

P 

I 

iThis  parameter  includes  the  combined  effect  of  DOB  and 

The  takeoff  angle  correction  is  small  for  teleseismic 
distances . 
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(a)  Direct  P Wave  (P) 


» 


(b)  Reflected  P Wave  (pP) 
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Figure  3.4. 


The  construction  of  the  first  arrival  P wave 
train  from  the  direct  and  reflected  phases. 
DOB  = 300  m,  = 2.4  km/sec,  tj  = 0.24  sec 


»l«t  II. tl 


.«  i.o  i.s  i.o  i.o 


t|»c 

(b)  t^  = 0.73  sec,  A = 18.8  my,  T = 1.4  sec. 

Figure  3.5.  P + pP  for  a series  of  delay  times  tL.  The  m^ 
measurement  using  standard  procedures  is  indi- 
cated. 


Figure  3. 


(c) 


(continued) 


= 0.97  sec,  A = 14.5  my,  T = 0.65  sec. 


= 1.2  sec,  A = 14.5  my,  T = 0.65  sec 


(d)  tr 


In  Table  3.1  the  data  of  Figures  3.4  and  3.5  are  written 
in  terms  of  standard  mb  « log  A/T  + B where  B is  the 
Gutenberg-Richter  distance  correction.  For  4087  km  a value 
of  3.53  was  used  for  B. 

Note  that  both  the  amplitude  A and  measurement 
period  T vary  with  lag  time.  Beyond  a certain  lag  time  the 
two  pulses  cease  to  interfere  and  mb  remains  constant  with 

depth. 


TABLE  3.1 


Figure 

Depth  (DOB) 
(meters) 

tL 

(sec) 

A 

(mu) 

T 

(sec) 

mb 

3.4 

300 

0.24 

20.7 

0.75 

4.97 

3.5(a) 

600 

0.49 

21.6 

0.95 

4.89 

3.5(b) 

900 

0.73 

18.8 

1.40 

4.66 

3.5(c) 

1200 

0.97 

14.5 

0.65 

4.88 

3.5(d) 

1500 

1.20 

14.5 

0,65 

4.88 

From  the  tabulation  (3.1)  we  observe  that  the  maxi- 
mum variation  in  mb  is  approximately  0.3  magnitude  units. 
This  variation  clearly  arises  from  the  procedure  employed 
in  measuring  mb<  That  is,  the  amplitude  and  effective 
period  used  are  chosen  to  be  those  associated  with  the 
largest  amplitude  cycle  within  the  first  two  or  three 
cycles  of  the  P wave  train.  For  very  shallow  events  the 
P and  pP  signals  arrive  with  relatively  small  time  separa- 
tion (compared  to  the  apparent  periods)  and  the  interference 
is  such  that  the  second  cycle  of  motion,  corresponding  to 
a sum  Of  P and  pP,  is  larger  than  the  first  cycle  (pure  P 


wave)  and  the  magnitude  measurement  is  therefore  increased, 
reflecting  a pP  contribution.  For  deeper  events  the  pP 
arrival  is  sufficiently  delayed  so  that  little  or  no  con- 
structive overlap  occurs,  and  the  magnitude  measurement 
is  then  made  from  the  first  cycle  of  motion.  In  this  latter 
case  the  magnitude  measurement  is  made  from  the  direct  P 
wave  alone,  while  in  the  former  case,  it  is  made  from  the 
constructively  summed  P and  pP  phases. 

If  the  magnitude  were  always  measured  from  the  first 

cycle  of  motion,  then  pP  interferrence , and  hence  depth  of 

burial,  would  have  little  or  no  effect  on  the  m,  measure- 

b 

ment.  (This  suggests  a change  in  the  mb  measurement:  that 

is  to  always  use  the  first  cycle  amplitude  of  the  P wave 
train  for  the  magnitude  measurement.  This  would  avoid  the 
pP  depth  controlled  contamination  and  also  any  sP  contamina- 
tion which  might  arise  when  tectonic  release  effects,  or 
other  S wave  generating  effects,  are  large.  The  latter  are 
discussed  in  parts  3.2  and  3.3  of  this  section.) 

There  are  three  ways  in  which  a change  in  source 
depth  affects  the  surface  wave  magnitude  M . Two  of  these 
have  already  been  discussed,  since  they  involve  changes  in 
near  source  material  properties  with  depth.  They  are,  first 
the  change  in  the  confining  pressure  with  depth  which  changes 
the  effective  material  strength  and  hence  ij>(<*>)  (see 
Section  II).  The  second  is  related  to  the  fact  that 
Ms  * log  [y  iK00)]  and  since  y,  the  elastic  ridigity,  will 
change  with  depth,  generally  increasing  with  increasing 
depth,  then  Mg  will  also  change.  Of  course,  other  material 
properties  will,  in  general,  change  with  depth  as  well,  but 
these  changes  would  be  reflected  in  changes  in  ipf00),  which  have 
been  discussed  in  Section  II  in  terms  of  the  rear  uurce 
material  properties. 
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Finally  the  M value  will  change  with  depth  due  to 
the  changes  in  the  depth  dependent  Rayleigh  wave  excitation 
factor.  However,  for  the  range  of  depths  involved  in 
underground  explosion  tests,  the  changes  in  Rayleigh  wave 
excitation  are  insignificant. 


3.3  SUMMARY  --  DEPTH  EFFECTS 


Since  m,  a log 


iK00) 

Tr 


and  Ms  log  [p  ^(°°)],  then 
U and  V 

P 

Further,  depth 


depth  changes  resuit  in  changes  in 
and  hence  variations  in  mb  and  Mg. 
changes  will  directly  change  P , the  confining 
pressure,  as  well  as  other  parameters 

controlling  <M°°)  variations  in  mb  and  Mg . 

The  relative  arrival  time  of  P,  pP  and  sP  vary 
with  depth.  For  explosions,  the  first  two  to 
three  cycles  of  the  P wave  train  usually  is  a 
composite  of  these  phases.  From  this  effect  mb 
as  usually  measured,  may  vary  by  roughly  0.3 
magnitude  units  for  source  depth  variations 
over  the  range  typical  for  underground  explo- 


sions . 


3.4  THE  EFFECT  OF  TECTONIC  RELEASE  ON  mb  AND  Mg 

The  theory  of  the  dynamic  relaxation  of  tectonic  stress 
due  to  the  creation  of  a cavity  and  surrounding  zone  of  frac- 
tured material  in  a prestressed  geologic  formation  was  deve- 
loped by  Archambeau  [1972]  to  explain  the  anomalous  shear 
waves  observed  from  underground  nuclear  explosions.  Tectonic 
release  adds  a double  couple  component  to  the  source  dis- 
placement field.  The  important  parameters  governing  the 
tectonic  release  contribution  to  the  source  are  the  amplitude 
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and  orientation  of  the  prestress,  the  radius  (R^)  of  the 
shatter  zone  created  by  the  explosion  shock  wave  and  the 
effective  size  (Rg)  of  the  prestressed  region  around  the 
explosion  cavity.  The  displacement  spectral  amplitudes 
scale  directly  with  prestress  and  with  the  cube  of  the 
shatter  zone  radius.  Since  R is  expected  to  scale  as 
W1/3,  the  latter  is  equivalent  to  tectonic  release  scaling 
with  yield.  The  parameter  Rg  controls  the  long  period 
portion  of  the  displacement  spectrum.  If  Rg  is  very 
large,  the  spectrum  is  essentially  flat  while  smaller  Rg 
causes  a falloff  in  the  long  period  level  and  a peaked 
spectrum. 

By  means  of  an  illustrative  example,  the  kinds  of  ef- 
fects tectonic  release  can  have  on  the  first  few  seconds 
of  the  P wave  train  will  now  be  demonstrated.  In  Figure 
3.6  two  examples  of  the  addition  of  tectonic  release  to  the 
record  of  Figure  3.5a  are  shown.  The  explosion  is  65  kT  in 
tuff  using  Source  T-3.  For  this  explosion  the  computed 
final  cavity  radius  was  54  meters.  The  shatter  zone  for 
the  tectonic  release  was  taken  to  be  'v  5-1/2  times  this 
radius  (300  m) . The  parameter  Rg  has  almost  no  effect 
on  P wave  ground  motion  as  measured  by  an  LRSM  short  period 
instrument.  The  convention  for  defining  the  prestress  and 
azimuth  are  indicated  in  Figure  3.7.  Note  that  the  addi- 
tion of  a vertical  component  to  the  prestress  field  has  a 
rather  large  effect  on  the  record. 

In  Figure  3.8  two  more  examples  of  tectonic  release 
are  shown.  In  these  examples  the  prepress  field  is  only 
slightly  perturbed  from  the  horizontal.  The  average  source- 
surface  travel  timos  for  both  P and  S waves  are  important 
when  tectonic  release  is  present.  The  effect  of  changing 
these  parameters  is  also  illustrated  in  Figure  3.8.  Of 
course,  these  two  figures  may  also  be  interpreted  to  show 
the  effect  of  changing  depth  of  burial. 
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The  effect  of  tectonic  release  on  the  first 
mantle  arrival  P wave  train  (vertical  compo- 
nent) . Explosion  source  T-3  (tu££) H*]? 

= 65  kT,  R0  = 300  m.  RS  - J H " £,  km 

VpS  = 2.4  km/sec,  V§S  = l km/sec,  A - 4087  km, 

azimuth  45°. 


MW 


iW 


C J 


(a)  Vfs 


r 


2.4  km/sec 
1.0  km/sec 


(b)  v®s 


1.6  km/sec 
0.8  km/ sec 


Figure  3.8.  The  effect  of  tectonic  release  on  the  first 
mantle  arrival  P wave  train  (vertical  compo 
nent).  All  explosion  and  tectonic  release 
parameters  are  identical  to  those  of  Figure 
3,7  except  S, , = 24  bars  for  both  (a)  and 
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While  in  none  of  these  examples  does  tectonic  release 
very  strongly  affect  the  measurement  of  m^,  it  is  important 
to  note  that  the  sP  phase  (which  may  be  large)  does  arrive 
within  the  first  few  seconds.  Thus,  it  may  strongly  in- 
fluence m,  measurements  under  certain  conditions,  dependent 

D 

on  the  magnitude  and  orientation  of  the  prestress  and  the 
net  time  delay  of  the  phase  with  respect  to  P. 

The  influence  that  tectonic  release  may  have  on  the 
Rayleigh  wave  due  to  an  explosion  is  demonstrated  by  comput- 
ing synthetic  seismograms  for  the  radial  component  Rayleigh 
wave  (measured  by  an  LRSM  instrument)  generated  by  several 
of  the  events  for  which  P-wave  seismograms  were  presented 
earlier. 

In  Figure  3.9  three  such  synthetic  seismograms  are 
shown.  Figure  3.9a  is  a Rayleigh  wave  at  2000  km  due  to  the 
65  kT  explosion  in  tuff  for  which  body  wave  seismograms  are 
shown  in  F gures  3.4  and  3.5.  In  Figures  3.6  and  3.8  the 
effect  of  two  slightly  different  tectonic  release  fields 
on  the  first  arrival  P wave  train  were  shown.  The  surface 
waves  for  these  fields  are  shown  in  Figures  3.9b  and  c. 

The  presence  of  this  level  of  tectonic  release  has  a 
dramatic  effect  on  the  amplitude  of  the  Rayleigh  wave  which 
is  increased  by  nearly  an  order  of  magnitude.  Note  that 
the  addition  of  a vertical  component  to  the  prestress  field 
has  little  effect  on  . Of  course,  the  contribution  of 
tectonic  release  to  the  Rayleigh  wave  is  azimuthally  depen- 
dent. This  dependence  is  shown  in  Figure  3.10  where  Ms 
enhancement  due  to  tectonic  release  characterized  by  the 
parameters  listed  in  Figure  3.8  is  plotted  versus  azimuth. 

Here  M enhancement  is:  M due  to  the  explosion  plus  tecto- 

s s . , 

nic  release  minus  Mg  due  to  the  explosion  alone.  In  such 

a plot  the  effect  of  the  source  receiver  travel  path  is 

removed.  It  is  seen  that  the  Mg  enhancement  is  0.7  to  1.0 

units  over  much  of  the  radiation  pattern. 
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(a)  Pure  explosion. 


(b)  Explosion  plus  tectonic  release:  Prestress  S12  = 120  bars 


Figure  3.9. 


Synthetic  seismograms  for  the  radial  component 
of  the  Rayleigh  waves  as  recorded  by  an  LRSM 
long  period  seismogram  at  an  epicentral  distance 
of  2000  km  and  an  azimuth  of  45°.  The  explosion 
was  modeled  by  the  tuff  source  T-3  at  a yield  of 
65  kT.  The  shatter  zone  radius  Rq  = 300  m O 5 
times  final  cavity  radius)  and  Rs  = 3 km. 


Figure  3.10.  The  azimuthal  distribution  of  Ms  enhancement 
due  to  tectonic  release  (specified  by  the 
parameters  of  Fig.  3.9c).  On  the  figure,  1.0 
represents  the  case  in  which  Ms  due  to  the  ex 
plosion  plus  tectonic  release  is  1.0  units 
greater  than  Ms  due  to  the  spherically  symmet 
ric  explosion  alone. 
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Finally,  it  should  be  mentioned  that  the  addition  of 
tectonic  release  can  account  for  the  large  Love  waves  often 
observed  from  underground  explosions.  Measurement  of  a 
surface  wave  magnitude  based  on  Love  waves  is  a direct 
measure  of  this  effect,  since  spherically  symmetric  explo- 
sions should  generate  no  Love  waves. 

3.5  SUMMARY  --  TECTONIC  RELEASE  EFFECTS 

Tectonic  release  effects  on  m^  and  Mg  include  the 
following : 

1.  A free  surface  reflected  sP  phase  is  generated 
which  appears  within  the  first  few  seconds  of 
the  teleseismic  P wave  record.  This  phase 
certainly  adds  complexity  to  the  record  but  al- 
so may  influence  mb  under  certain  circumstances, 
in  particular  when  the  medium  is  highly  stressed, 
has  a high  shear  velocity  and/or  the  source  is 
shallow. 

2.  Love  waves  are  generated. 

3.  The  Rayleigh  wave  may  be  greatly  increased  from 
its  value  for  a similar  explosion  with  no  tecto- 
nic release. 

4.  All  tectonic  release  effects  are  strongly  depen- 
dent on  the  orientation  and  magnitude  of  the 
prestress  field,  and  hence  are  strongly  azimuthally 
dependent . 

3.6  THE  EFFECT  OF  SOURCE  ASYMMETRY  ON  Mg  AND  mb 

Various  factors  characterizing  the  explosion  may  be 
responsible  for  adding  a double  couple  (quadrupole)  pertur- 
bation to  the  basic  explosion  center  of  dilatation  (monopole). 
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These  may  include  jointing  and  other  source  region  asymme- 
tries, tectonic  stress  release  due  to  the  introduction  of 
a large  cavity  in  a prestressed  media,  earthquake  trigger- 
ing and  the  deliberate  introduction  of  asymmetry  via  emplace- 
ment cavity  design.  The  latter  possibility  has  been  dis- 
cussed by  Cherry,  et  al.  [1974]  in  a theoretical  and  ex- 
perimental study. 

Using  a computational  technique  in  which  Eulerian 
hydrodynamic  codes  were  linked  to  a large  displacement,  non- 
linear Lagrangian  code  and  finally  to  elastic  wave  propaga 
tion  codes,  Cherry,  et  al.,  computed  the  teleseismic  ground 
motion  generated  by  a nuclear  explosion  in  a tunnel  of  as- 
pect ratio  200  x 3 in  tuff.  The  resulting  free  field  dis- 
placement radiation  pattern  as  computed  10  km  from  the 
source  (retaining  only  far  field  terms)  is  shown  in  Figure 
3.11,  The  strong  double  couple  perturbation  of  the  dis- 
placement field  is  evident,  especially  in  the  20  second 
pattern. 

From  the  displacement  field  pattern  of  Figure  3.11, 
Cherry,  et  al.,  computed  surface  waves  at  teleseismic  dis- 
tances. The  effect  of  the  layered  earth  structure  was  re- 
moved by  normalizing  the  resulting  surface  wave  amplitudes 
to  those  resulting  from  the  center  of  dilatation  component 
of  the  source.  For  the  tunnel  oriented  at  45°  to  the 
horizontal,  the  resulting  enhancement  of  Mg  is  shown  in 
Figure  3.12.  Averaged  azimuthally,  Mg  due  to  this  emplace- 
ment configuration  was  found  to  be  0.8  units  larger  then 
that  due  to  the  center  of  dilatation  component  alone.  On 
the  other  hand,  for  the  tunnel  oriented  vertically,  the 
M enhancement  was  0.3  units  while  for  the  horizontal  orienta 
tion  the  average  M enhancement  was  -0.4  units  (Mg  was  less 
than  that  due  to  the  center  of  dilatation  component).  For 
all  orientations  body  wave  magnitudes  was  essentially  the 
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(a)  1 Second  radiation  pattern 


(b)  20  Second  radiation  pattern 


Figure  3.11.  Far  field  displacement  radiation  pattern  10  km 
from  the  hypocenter.  Both  1 and  20  second  patterns  are 
shown  and  the  ratio  of  S to  P wave  displacements  are  indicated 
on  the  right  half  of  each  pattern.  The  P wave  radiation  from 
the  monopole  alone  is  also  indicated. 


same  as  that  from  the  spherically  symmetric  portion  of  the 
source.  This  is  apparent  from  the  illustration  of  the  1 
.,\r  n<  ’iation  pattern  in  Figure  ^.11. 

While  further  verification  (the  extremely  complex 
calculation  was  carried  out  only  once)  of  the  quantitative 
results  presented  above  is  required  before  complete  con- 
fidence is  warranted,  qualitatively  it  seems  clear  that  a 
peiturbation  of  the  source  from  a pure  center  of  dilatation 
may  strongly  influence  surface  wave  amplitudes.  Further, 
this  effect  is  quite  dependent  on  the  orientation  of  the 
source  asymmetry. 

Another  method  by  which  significant  source  asymmetries 
may  be  introduced  is  the  detonation  of  two  or  more  explosions 
simultaneously  or  in  sequence.  The  shock  interaction  of  two 
nearby  explosions  detonated  simultaneously  is  expected  to 
generate  elastic  waves  similar  to  those  from  an  explosion 
in  a tunnel  (discussed  above).  On  the  other  hand,  a suite 
of  explosions  may  be  generated  in  sequence  in  such  a way  as 
to  simulate  many  features  of  an  earthquake.  This  constitutes 
an  evasion  procedure  which  has  received  some  study  (e.g., 
Archambeau,  ACDA  report,  California  Institute  of  Technology, 
1974)  and  it  appears  that  the  mb  and  M values  of  the  sum 
o che  events  of  the  series  can  be  determined  by  use  of 
ultra  narrow  band  filtering  of  both  the  P wave  and  surface 
wave  trains  from  the  multiple  event.  Because  the  relation- 
ship of  mb  to  Mg  for  explosions  is  different  ‘han  for  earth- 
quakes, the  event  appears  as  an  explosion  under  an  mb-Ms 
discrimination  criteria.  In  terms  of  a yield  estimate  for 
such  a multiple  event  using  mb  and  Mg;  one  would  obtain  an 
estimate  of  the  total  yield  of  all  the  explosions  of  the 
series . 

Finally,  s ich  source  environment  asymmetries  as 
jointing  may  be  responsible  for  many  of  the  anomalous 
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radiation  pattern,  Love  wave  geneiation  and  shear  wave 
generation  effects  observed  from  underground  nuclear  explo- 
-i.  ...  ' .’ith  other  asymmetric  Ai 

this,  in  general,  is  not  expected  to  have  much  effect  on 
mb  but  may  have  large  effects  on  Rayleigh  wave  and  Love 
wave  amplitudes.  Hence,  yield  estimates  derived  from  mfa 
would  be  relatively  unaffected. 


3 . 7 SUMMARY 


1.  Source  asymmetry  due  to  emplacement  ca\ity  design 
produces  effects  similar  to  tectonic  release  with 
with  the  perturbation  of  mb  and  Mg  being  of  the 
order  of  that  due  to  tectonic  release  with  modest 
(75  bars)  prestress.  Mg  levels  car.  fluctuate  by 
as  much  as  a magnitude  unit  depending  on  cavity 
orientation.  On  the  other  hand,  mb  variation 

from  that  due  to  a symmetric  source  is  much  smaller, 
on  the  order  of  0.1  magnitude  units. 

2.  Jointing  and  other  medium  anisotropies  are  thought 
to  have  effects  on  mb  and  Ms  which  may  be  of  the 
same  order  of  magnitude  as  those  produced  by 
emplacement  cavity  design. 

3.  Multiple  explosions  apparently  can  be  discriminated 
using  revised  mb  and  Mg  measures.  The  yield  esti- 
mate obtained  from  these  revised  magnitude  values 
would  be  (roughly)  equal  to  the  sum  of  the  yields 
from  all  the  events  of  the  series. 
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IV.  DEPENDENCE  OF  MAGNITUDE  MEASUREMENTS  ON 
VARIATIONS  IN  EARTH  STRUCTURE 

Reliable  predictions  of  the  teleseismic  signature 
of  body  and  surface  waves  excited  by  an  underground  explosion 
depend  strongly  on  the  knowledge  of  a valid  velocity  and  Q 
model  for  the  source-receiver  earth  structure.  Variations 
in  the  velocity  and  Q structure  of  the  earth's  crust  and 
upper  mantle  give  rise  to  a complex  seismogram.  These 
variations  must  be  properly  modeled  before  statistically 
reliable  determinations  of  body  wave  magnitudes  (mb)  and 
surface  wave  magnitudes  (Mg)  can  be  obtained.  In  this  section 
of  the  report  we  will  emphasize  the  influence  of  the  inter- 
vening earth  structure  on  both  the  body  wave  and  surface 
wave  signatures  and  relate  these  effects  to  magnitude  deter- 
minations . 

4.1  BODY  WAVES 

The  standard  Gutenberg- Richter  formula  for  body  wave 
magnitude  may  be  written 


mb  = log [A/T]  + B , 


(4.1) 


where 


Maximum  zero- to-peak  amplitude  (millimicrons) 
of  the  most  prominent  P phase  in  the  first 
few  cycles, 

Period  of  cycle  on  which  A is  measured, 

Station  correction;  as  usually  applied  (e.g., 
by  the  LRSM  network)  B corrects  for  distance 
only,  implicitly  assuming  a laterally  homo- 
geneous earth. 
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Perhaps  the  most  prominent  characteristic  of  compila- 
* tions  of  observed  mb  data  based  on  Eq.  (4.1)  is  its  wide 

scatter  por  instance,  the  standard  deviation  associated 
with  LRSM  mb  values  for  U.  S.  underground  explosions  is 
typically  as  much  as  0.4  or  0.5.  As  an  example,  the  LRSM 
o data  for  the  Rulison  event  taken  from  the  Teledyne  Geotech 

shot  report  is  shown  in  Figure  4.1.  The  magnitude  assigned 
to  this  event  was  based  on  data  from  the  eleven  stations  at 
A > 17°  and  was  mb  - 4.62  ♦ 0. 36.  The  individual  station 
magnitudes  varied  from  4.07  to  5.21,  a range  of  more  than  a 
magnitude  unit. 

The  Rulison  explosion  was  detonated  in  western  Colorado. 
However,  the  state  of  the  mb  data  is  no  different  for  NTS 
events.  For  example,  the  LRSM  data  for  Piledriver  is  shown 
in  Figure  4.2.  In  this  case  the  average  mb  is  5.43  + 0.42 
for  A > 16°  and  individual  values  for  the  ten  stations  in  this 
region  vary  from  4.69  to  6.05,  a range  of  1.36  units. 

There  is  no  reason  to  believe  that  any  other  seismic 
network  offers  a less  scattered  mb  data  set  than  the  LRSM. 
Possible  reasons  for  this  variation  that  may  be  cited  are 
many  and  include  such  correctable  causes  as  calibration 
and  measurement  inconsistencies.  However,  by  far  the  most 
important  are  path  effects.  That  is,  the  amplitude  and 
even  the  frequency  content  of  the  first  motion  P wave  at 
individual  stations  of  any  seismic  network  will  exhibit 
wide  variation  due  to  fundamental  properties  of  the  earth. 

The  P wave  from  an  explosion  which  reaches  several 
teleseismic  stations  has  traveled  a different  path  for  each, 
less  so  in  the  lower  mantle,  more  so  near  the  surface.  The 
effects  of  local  and  regional  receiver  geology  on  the 
characteristics  of  signals  recorded  are  depicted  in  Figures 
4.3a  and  4.3b.  Figure  4.3a  shows  the  locations  of  five  portable 
LRSM  stations  that  were  located  within  4 km  of  the  same 
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Figure  4.3(a).  Locations  of  LRSM  portable  systems  for  the 

Corduroy  Event.  (AFTAC/Teledyne  Shot  Report, 
1966.) 
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Figure  4. 3(b) . 


p phases  from  the  Corduroy  Event  recorded 
by  the  LRSM  portable  systems,  normalized  to 
the  same  magnification  at  each  site.  All 
traces  are  short-period  vertical  at  a mag- 
nification of  56.2  (at  1 cps)  . (AFTAC/Telc- 
dyne  Shot  Report,  1966.) 
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distance  from  the  shot  site.  While  all  the  recording  sites 
are  located  in  the  eastern  United  States,  each  is  situated 
in  a region  of  distinctly  different  geologic  character 
ranging  from  Pre-Cambrian  basement  granite  at  AX-AL  to  a 
section  of  unconsolidated  alluvium  and  glacial  drift  at  ED-MI. 
Figure  4.3b  shows  the  short-period  vertical  component  of  the 
P wave  signals  recorded  at  each  of  these  five  locations  from 
the  Corduroy  explosion  at  NTS.  All  traces  are  at  the  same 
magnification.  The  most  significant  result  is  that  mb  (AFTAC/ 
Teledyne  Shot  Report,  1966)  ranges  from  5.2  at  RS-KY  to  5.7 
at  ED-MI  and  B0-A1,  a spread  of  0.5  mb  units. 

In  addition  to  the  near  station  structure  we  must  also 
consider  the  source  to  station  travel  path.  Body  waves  propa- 
gating through  a radially  stratified  earth  structure  undergo 
radial  spreading,  attenuation  and  dispersion.  Using  ray 
theory,  a complex  transfer  function  that  accounts  for  these 
effects  can  be  written  as  follows; 

Ap(w)  - Ap  exp  ^ [1  * -jr  *n  ] (4.2) 


where 


A 


P 


the  multiplication  factor  due  to  radial 
spreading , 


Qp  = the  dissipation  constant  for  the  ray  path, 

= the  high  frequency  cutoff  parameter  (Strick, 
1970). 


The  amplitude  of  the  transfer  function  |Ap(w)|  at 
u = 2ir  (T  = 1 sec)  for  the  first  P arrival  is  plotted  versus 
distance  in  Figure  4.4  for  three  different  earth  models.  Each 
transfer  function  was  calculated  for  a source  depth  of  1.2  km. 
For  a particular  earth  model  the  transfer  functions  for  other 
frequency  components  of  interest  for  teleseismic  P waves 


i 
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g exhibit  a similar  distance  dependence.  Also,  since  the 

amplitude  of  the  transfer  function  |Ap(2n)2  is  plotted 
along  trie  ordinate  on  a logarithmic  ^tun.  c..c  aic 

proportional  to  log  A/T  in  Eq . (1)  for  T = 1 second.  Thus 
the  linear  scale  along  the  ordinate  in  Figure  4.4  gives 
increments  in  terms  of  m^. 

The  dashed  curves  in  Figure  4.4  give  the  transfer 
functions  for  two  contrasting  earth  structures.  The  curve 
w labeled  CIT  109A/P  0.5Q  corresponds  to  an  earth  structure 

appropriate,  in  terms  of  the  upper  mantle  and  crustal 
structure  in  the  upper  200  km,  to  the  Basin  and  Range 
province  of  the  western  United  States  (Archambeau,  e t al . , 1969) . 
G The  curve  labeled  Masse'  (74)  High  Q corresponds  to  an  earth 

structure  more  appropriate  to  the  crust  and  upper  mantle 
beneath  central  and  eastern  North  America.  In  general 
these  two  models  correspond  to  the  crust-upper  mantle 
structure  in  a tectonic  province  (CIT  109  A/P  0.5Q)  and  a 
shield  province  (Masse'  (74)  High  Q) . 

In  terms  of  m^  the  mean  separation  of  these  two  trans- 
fer functions  is  approximately  0.6  units  for  distances  greater 
than  2500  km.  In  the  distance  range  1600  to  '.500  km  the 
separation  increases  with  decreasing  distance  from  0.6  units 
at  2500  km  to  1.7  units  at  1600  km.  This  latter  distance 
range  is  where  the  major  differences  in  the  two  models, 
namely  the  top  200  km  of  structure,  arc  most  pronounced. 
Considering  each  transfer  function  separately,  the  mean 
varies  by  0.2  to  0.5  m}  units  between  2500  and  7000  km, 
with  a larger  variation  of  0.4  m^  units  occurring  in  the 
distance  range  near  4000  km  in  both  models. 

So  far,  we  have  treated  the  source-receiver  earth 
structure  as  being  laterally  homogeneous.  Ihc  curve  labeled 
I "Split  Structure"  in  Figure  4.4  is  an  attempt  to  approximate 

a more  realistic  earth  model  for  a particular  experiment 
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(NTS  shots  recorded  at  LRSM  stations  in  eastern  North  America). 
Specifically,  the  split  structure  is  appropriate  to  an  earth 
nt-Jelc  two  velocity/attenuation  (Q)  models;  one  for  the 

region  from  t lie  source  to  the  ray's  deepest  penetration  (CIl 
109  A/P,  0.5Q)  and  one  for  the  region  from  the  point  of  deepest 
penetration  to  the  receiver  (Masse'  (74)  High  Q) . Note  that 
at  distances  greater  than  3000  km  the  mean  of  the  split  struc- 
ture transfer  function  varies  by  0.2  to  0.3  m^  units  while 
at  shorter  distances  it  reflects  the  contrast  in  the  two 
basic  eartli  structure  models. 

The  importance  of  the  earth's  structure  on  m^  can  be 
demonstrated  in  another  way.  In  Figure  4.5  log  |Ap(2n)| 

(given  by  the  solid  curve  and  log  A/T  = m^  -B  from  Ilq.  (4.1) 
with  T = 1 sec  (given  by  the  dashed  curve)  are  plotted  versus 
distance.  A value  of  m^  = 4.4  was  chosen  in  order  to  scale 
the  dashed  curve. 

As  previously  mentioned,  B corrects  foi  distance  only, 
assuming  a laterally  homogeneous  earth.  While  Jie  two  curves 
are  quite  close  over  the  distance  range  5000  to  5000  km,  they 
diverge  by  0.5  m^  units  or  more  on  either  end  of  this  dis- 
tance range.  Clearly,  to  reduce  the  variation  in  m^  with 
distance,  B must  be  dependent  on  the  particular  source  region 
from  which  the  signal  emanates,  the  source-receiver  earth 
structure  and  the  regional  structure  in  the  vicinity  of  each 
receiver  in  a particular  network. 

By  a careful  determination  of  station  corrections, 
it  is  possible  to  significantly  reduce  the  standard  devia- 
tion and  range  of  an  m^  data  set.  Such  an  exercise  has  been 
carried  out  by  Basham  (1969)  on  the  Canadian  data  from  N S 
explosions.  He  was  able  to  reduce  the  average  standard 
deviation  for  m^  due  to  eleven  explosions  (with  mean  m^ 
ranging  from  4.46  to  6.29)  to  0.20.  With  the  passage  of 
time  and  accumulation  of  data,  it  is  reasonable  to  expect 
further  improvements. 


MHMMp 


The  important  point  to  be  borne  in  mind  regarding  m^ 
is  that  it  is  not  only  method  dependent  (i.e.,  which  peak 
is  meas.ieu,  should  che  amplitude  be  div'd*'  ’ \:J  t 'tc.) 
but  also  dependent  on  the  stations  themselves;  tne  character- 
istics of  the  sites  on  which  they  are  located  and  the  source- 
station  travel  path.  Much  statistical  analysis  of  data,  like 
that  done  by  Basham  (1969)  remains  to  be  done  before  reliable 
"mean"  or  "average"  m^  values  can  be  given  for  past  or  future 
events . 

4.2  SUMMARY 

1.  m^  variations  due  to  earth  structure,  including 
variations  of  curst-upper  mantle  structure  at 
different  source  sites,  are  typically  of  the  order 
of  0.5  magnitude  units.  Examples  of  variations 
greater  than  1 magnitude  unit  are  not  uncommon 
however.  Current  models  of  the  earth,  plus  tak- 
ing account  of  local  station  structure  and  source 
site  structure  can  probably  reduce  the  station 
to  station  variation,  after  correction,  to  from 
0.1  to  0.2  magnitude  units  and  provide  an  average 
m^  which  is  accurate  to  at  least  0.1  magnitude 
unit  and  probably  much  better  than  this  if  only 
random  station  to  station  fluctuations  remain 
after  correction. 

4.5  SURFACE  WAVES 

The  amplitude  and  frequency  of  surface  waves  recorded 
at  teleseismic  distances  from  underground  explosions  are 
quite  sensitive  to  variations  in  the  source-receiver  velocity 
(P,S)  profile  of  the  crust  and  upper  mantle.  The  effect  of 
such  variations  may  be  seen  by  the  following  examples. 
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Synthetic  seismograms  at  4000  km  from  the  source  were  gene- 
rated for  three  earth  structure  models.  The  techniques 
dwscrik  ' ~y  Cherry,  e_t  (1974)  were  '•*=;- 1 ’ ” * v ~ calr-’a 

tions.  The  three  models,  which  are  presented  in  Figure  4.0, 
are  representative  of  the  earth  sturcture  in  the  vicinity  of 
NTS  and  are  therefore  quite  similar.  The  synthetic  seismo- 
grams are  appropriate  to  the  Piledriver  event  which  was  a 
61  kt  explosion  at  a burial  depth  of  463  meters  in  NTS 
granite . 

The  earth  structure  affects  both  the  arrival  time 
(group  velocity)  and  amplitudes  of  the  various  frequency  com- 
ponents. Surface  waves  from  different  events  are  often  com- 
pared through  the  surface  wave  magnitude,  M . The  form  of 

M used  here  is 
s 

M » log  (Jii)  ♦ 1.656  logA  + 3.119  (4.3) 

where  A]f  is  the  maximum  zero-to-peak  amplitude  in  microns 
of  the  horizontal  component  of  the  Rayleigh  wave  having  a 
period  T near  20  seconds,  and  A is  the  epicentral  dis- 
tance in  degrees.  This  formula  is  to  be  applied  only  in  the 
distance  range  15°  < A < 130°. 

In  Table  4.1  the  M values  based  on  the  synthetic 

s 

seismograms  are  summarized.  The  seismograms  were  made  to 
correspond  to  those  written  by  an  LRSM  long  period  seismo- 
meter. From  the  table  we  see  that  the  earth  structure  can 
have  a noticeable  effect  on  ?1  , even  for  these  three  similar 
models.  A more  significant  variation  in  M results  from  the 
choice  of  the  phase  of  maximum  amplitude,  the  Airy  phase, 

for  M measurements,  rather  than  the  amplitude  at  20  seconds, 
s 

Values  of  Mg  from  the  two  are  really  not  comparable. 

The  extent  of  the  Ms  variations  in  Table  4.1,  namely 
0.2  units,  should  be  considered  a minimum  in  view  of  the 
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similarity  of  the  three  cartli  structure  models  (Figure  4.0). 
Propagation  of  surface  waves  through  cont inental -oceanic 
V"  ndr ~ : *■'  ">r  into  regions  oF  quite  different  crustal  struc- 
ture can,  of  course,  result  in  average  velocity-depth  pro- 
files quite  different  from  those  of  Figure  4.0  and  introduce 
a more  significant  variation  into  the  M measurements. 

An  example  of  the  wide  scatter  in  M determinations 
associated  with  a source  located  in  one  region  and  the  re- 
cording stations  in  different  regions  (in  terms  of  velocity- 
depth  profiles)  is  given  in  Figure  4.7.  In  particular,  we 
see  that  the  M.  values  for  the  Rulison  explosion,  determined 
at  A > 15°,  scatter  by  as  much  as  1.4  units.  The  M measure- 
ments were  made  on  Rayleigh  waves  in  the  period  range  11  to 
IP  seconds  with  a preponderance  of  measurements  at  the  shorter 
periods.  The  average  M in  this  distance  range  is  given 
(Lambert  and  Ahner,  1972)  as  = 3.99  + 0.44  based  on  12 
stations  (Figure  4.7).  Since  a large  variation  in  M values 
(0.94  uni  is)  is  observed  at  stations  along  nearly  the  same 
azimuth  from  Rulison  and  oince  the  observed  Love  waves  were 
considerably  smaller  (1/2  to  1/3)  than  Rayleigh  waves  of 
similar  period,  the  wide  scatter  in  M is  most  probably  a 
result  of  path  effects  rather  than  possible  tectonic  strain 
release  or  source  asymmetries. 

Marshall  and  Basham  (1972)  argue  that  if  measurements 
of  M are  confined  to  Rayleigh  waves  of  20  second  periods, 
then  for  distances  > 25°  the  largest  relative  difference 
due  to  path  alone  should  be  about  0.2,  M units.  However, 
for  A < 25°,  it  is  necessary  to  make  careful  path  corrections 
for  each  source-receiver  combination.  That  is,  the  accuracy 
of  the  earth  structure  model  will  determine  the  reliability 
of  the  Ms  determination. 
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4 . 4 SUMMARY 


The  variations  in  observed  values  arc  to  some 
considerable  extent  due  to  measurements  ot  Ray- 
leigh waves  of  different  apparent  periods  at 
different  points  of  observations.  This  can  be 
avoided  by  measuring  only  the  20  second  Rayleigh 
wave  amplitude.  True  structural  variations  can 
nevertheless  cause  variations  in  point  to  point 
measurements  of  at  least  0.2  magnitude  units, 
and  in  some  cases  of  around  0.4  units.  Correc- 
tions for  structure,  amounting  to-  different  B 
value  corrections  for  different  paths  based  on 
currently  obtainable  crust-mantle  structure, 
should  reduce  this  variation  to  a (random)  + 0.1 
to  0.2  variation  and  a resulting  mean  Mg  value 
which  should  be  accurate  to  less  than  0.1  magni- 
tude unit. 
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V.  CONCLUSIONS 


Results  from  a deterministic  model  have  been  presented 
which  permit  a systematic  analysis  of  the  factors  which  affect 
body  wave  and  surface  wave  magnitudes,  m^  and  M , obtained 
from  a nuclear  explosion.  The  model  merges  nonlinear  shock 
code  calculations  in  the  near-field  with  seismic  code  calcu- 
lations of  far-field  stress  wave  propagation  through  a 
realistic  earth  structure. 

In  this  study  the  model  was  used  to  produce  theoretical 
seismograms  at  teleseismic  distances.  The  theoretical  seismo- 
grams were  then  used  to  determine  m^  and  Mg.  Variations 
in  these  magnitudes  have  been  correlated  with  variations  in  the 
parameters  in  the  model. 

The  effect  of  near  source  material  properties  on 
magnitude  were  presented  in  Section  II.  Magnitude  dependence 
on  depth  of  burial,  medium  prestress  and  source  asymmetry 
are  discussed  in  Section  III.  Effects  of  earth  structure 
on  teleseismic  magnitudes  are  given  in  Section  IV. 

The  important  conclusions  concerning  teleseismic 
magnitude  variations  are  contained  in  each  section  and  will 
not  be  repeated  here.  However,  we  do  wish  to  emphasize  that 
a complete  predictive  capability  now  exists  for  estimating 
body  wave  and  surface  wave  ground  motion  caused  by  an  explosive 
source  as  a function  of  depth  of  burial,  emplacement  technique, 
near-source  environment,  tectonic  stress  release  and  earth 
structure.  These  ground  motion  predictions  permit  an 
"inversion"  of  field  observations  in  order  to  uncover  device 
yield  and  also  provide  an  evaluation  of  the  effectiveness  of 
existing  and  proposed  discriminatory  criteria. 

Refinement  of  the  model  will  certainly  occur  in  the 
future.  These  modifications  will  only  proceed  as  a result 


® of  comparing  predicted  and  recorded  ground  motion  over  as 

large  a frequency  and/or  time  interval  as  possible.  This 
comparative  study  is  proceeding  at  S3.  The  current  status 
of  the  model  reflects  our  satisfaction  with  its  ability  to 
® match  the  data  available  to  us  at  the  present  time. 
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